Singularities and Poincar\'{e} Indexes of Electromagnetic Multipoles

TL;DR

This paper systematically maps the singularities and topological indexes of electromagnetic multipoles, revealing their connection to bound states in the continuum and unifying multipole theory with topology.

Contribution

It provides a comprehensive analysis of multipolar singularities and their indexes, linking them to topological charges and the formation of bound states in the continuum.

Findings

Index sum over the sphere is always 2, consistent with Poincaré-Hopf theorem.

Bound states in the continuum are related to overlapping multipolar singularities.

Multipolar singularities are connected to topological charges of bound states.

Abstract

Electromagnetic multipoles have been broadly adopted as a fundamental language throughout photonics, of which general features such as radiation patterns and polarization distributions are generically known, while their singularities and topological properties have mostly slipped into oblivion. Here we map all the singularities of multipolar radiations of different orders, identify their indexes, and show explicitly the index sum over the entire momentum sphere is always $2$, consistent with the Poincar\'{e}-Hopf theorem. Upon those revealed properties, we can attribute the formation of bound states in the continuum to overlapping of multipolar singularities with open radiation channels. This insight unveils a subtle connection between indexes of multipolar singularities and topological charges of those bound states, revealing that essentially they are the same. Our work has fused two fundamental and sweeping concepts of multipoles and topologies, which can potentially bring unforeseen opportunities for many multipole related fields within and beyond photonics.

Tables

Table 1. Table 1: Distributions of the singularities and their indexes for radiated vector fields from a multipole of order ( l , m ) 𝑙 𝑚 (l,m) .

Positions	Index	Number	Total Indexes
$\left|\cos \left(\theta \right)\right|=1$	$1-\left|m\right|$	2	$2\left(1-\left|m\right|\right)$
$\cos \left(m\varphi \right)=P_l^m\left(\cos \theta \right)$		$2\left|m\right|\times$	$-2\left|m\right|\times$
$=0,\left|\cos \left(\theta \right)\right|\ne 1$	-1 $\left(m\ne 0\right)$	$\left(l-\left|m\right|\right)$	$\left(l-\left|m\right|\right)$
$\sin \left(m\varphi \right)=\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }$	0 $\left(m=0\right)$	$2\left|m\right|(l-$	$2\left|m\right|\times$
$=0,\left|\cos \left(\theta \right)\right|\ne 1$	+1 $\left(m\ne 0\right)$	$\left|m\right|+1)$	$\left(l-\left|m\right|+1\right)$

Table 2. Table S1: Distributions of the singularities and their indexes for 𝕄 l ​ m e ​ ( | m | ≠ 1 ) subscript superscript 𝕄 𝑒 𝑙 𝑚 𝑚 1 \mathbb{M}^{e}_{lm}(|m|\neq 1) .

Positions	Index	Number	Total Indexes
$\left|\cos \left(\theta \right)\right|=1$	$1+\left|m\right|$	2	$2\left(1+\left|m\right|\right)$
$\cos \left(m\varphi \right)=P_l^m\left(\cos \theta \right)$		$2\left|m\right|\times$	$2\left|m\right|\times$
$=0,\left|\cos \left(\theta \right)\right|\ne 1$	+1 $\left(m\ne 0\right)$	$\left(l-\left|m\right|\right)$	$\left(l-\left|m\right|\right)$
$\sin \left(m\varphi \right)=\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }$	0 $\left(m=0\right)$	$2\left|m\right|(l-$	$-2\left|m\right|\times$
$=0,\left|\cos \left(\theta \right)\right|\ne 1$	-1 $\left(m\ne 0\right)$	$\left|m\right|+1)$	$\left(l-\left|m\right|+1\right)$

Equations

$${\mathbf{M}}_{lm}=-\frac{\mathrm{m}}{\sin \theta }\sin (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }$$

$$-\cos (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\varphi }.$$

$$\mathbf{Ind}=1-|m|,|\cos \theta |=1.$$

$$\mathbf{Ind}=-1,\cos (m\varphi )=P_l^m\left(\cos \theta \right)=0;$$

$$\mathbf{Ind}=+1,\sin (m\varphi )=\frac{d{P}_l^m(\cos \theta )}{d\theta }=0.$$

$${\mathbb{M}}_{lm}=\frac{\mathrm{m}}{\sin \theta }\sin (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }$$

$$-\cos (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\varphi }.$$

$$\mathbf{M}=▽\times (\mathbf{r}\Psi ),\mathbf{N}=\frac{▽\times \mathbf{M}}{k},$$

$$\begin{array}{r}{\Psi }_{lm}^e=\cos (m\varphi )P_l^m\left(\cos \theta \right)z_n(kr),\\ {\Psi }_{lm}^o=\sin (m\varphi )P_l^m\left(\cos \theta \right)z_n(kr),\end{array}$$

$$\begin{array}{r}\begin{array}{r}{\mathbf{M}}_{lm}^e=-\frac{\mathrm{m}}{\sin \theta }\sin (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }\\ -\cos (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\varphi },\end{array}\\ \begin{array}{r}{\mathbf{M}}_{lm}^o=\frac{\mathrm{m}}{\sin \theta }\cos (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }\\ -\sin (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\varphi },\end{array}\\ \begin{array}{r}{\mathbf{N}}_{lm}^e=-\frac{\mathrm{m}}{\sin \theta }\sin (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\varphi }\\ +\cos (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\theta },\end{array}\\ \begin{array}{r}{\mathbf{N}}_{lm}^o=\frac{\mathrm{m}}{\sin \theta }\cos (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\varphi }\\ +\sin (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\theta }.\end{array}\end{array}$$

$$|{\mathbf{M}}_{lm}^e|=0.$$

$$\mathbf{Ind}({\theta }_{\mathrm{s}},{\varphi }_{\mathrm{s}})=\frac{1}{2\pi }{\oint }_{\gamma }d\alpha ,$$

$${\oint }_{\gamma (\theta \to 0,\pi )}d\alpha =2\pi +2\pi |m|\sigma ,$$

$$\mathbf{Ind}(\theta =0,\pi )=1-|m|,$$

$${\mathbb{M}}_{lm}^e=\frac{\mathrm{m}}{\sin \theta }\sin (m\varphi )P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }$$

$$-\cos (m\varphi )\frac{d{P}_l^m\left(\cos \theta \right)}{d\theta }{\hat{\mathbf{e}}}_{\varphi }.$$

$$\mathbb{IND}(\theta =0,\pi )=1+|m|.$$

$$\cos (m\varphi )=P_l^m\left(\cos \theta \right)=0,$$

$$\sin (m\varphi )=\frac{d{P}_l^m(\cos \theta )}{d\theta }=0.$$

$$\mathbf{Ind}=-1,\cos (m\varphi )=P_l^m\left(\cos \theta \right)=0;$$

$$\mathbf{Ind}=+1,\sin (m\varphi )=\frac{d{P}_l^m(\cos \theta )}{d\theta }=0.$$

$$\underset{\mathbf{Poles}}{\underbrace{2(1-\left|m\right|)}}\underset{Eq.(S13)}{\underbrace{+2m(l-\left|m\right|+1)}}\underset{Eq.(S12)}{\underbrace{-2m(l-\left|m\right|)}}$$

$$\mathbb{IND}=+1,\cos (m\varphi )=P_l^m\left(\cos \theta \right)=0;$$

$$\mathbb{IND}=-1,\sin (m\varphi )=\frac{d{P}_l^m(\cos \theta )}{d\theta }=0.$$

$$\underset{\mathbf{Poles}}{\underbrace{2(1+\left|m\right|)}}\underset{Eq.(S16)}{\underbrace{-2m(l-\left|m\right|+1)}}\underset{Eq.(S15)}{\underbrace{+2m(l-\left|m\right|)}},$$

$$a_{lm}=\frac{{\left(-i\right)}^{l-1}{k}^2\eta Q_{lm}}{{\left[\pi \left(2l+1\right)\right]}^{\frac{1}{2}}}\underset{\mathbf{V}}{\int \int \int }\mathbf{J}\left(\mathbf{r}\right)\cdot {\mathbf{S}}_{lm}{d}^{3}r;$$

$$b_{lm}=\frac{{\left(-i\right)}^{l+1}{k}^2\eta Q_{lm}}{{\left[\pi \left(2l+1\right)\right]}^{\frac{1}{2}}}\underset{\mathbf{V}}{\int \int \int }\mathbf{J}\left(\mathbf{r}\right)\cdot {\mathbf{T}}_{lm}{d}^{3}r,$$

$$\mathbf{J}\left(\mathbf{r}\right)=-i\omega {\epsilon }_0\left[{\epsilon }_r\left(\mathbf{r}\right)-{\epsilon }_b\right]\mathbf{E}\left(\mathbf{r}\right),$$

$$\begin{array}{r}{\mathbf{S}}_{lm}=\exp \left(-im\varphi \right)\left[{\Pi }_l\left(kr\right)+{{\Pi }_l}^{\prime \prime }\left(kr\right)\right]P_l^m\left(\cos \theta \right){\hat{\mathbf{e}}}_r\\ +\exp \left(-im\varphi \right)\frac{{{\Pi }_l}^{\prime }\left(kr\right)}{kr}\left[{\tau }_{lm}\left(\theta \right){\hat{\mathbf{e}}}_{\theta }-i{\pi }_{lm}\left(\theta \right){\hat{\mathbf{e}}}_{\varphi }\right]\end{array};$$

$${\mathbf{T}}_{lm}=\exp \left(-im\varphi \right)j_l\left(kr\right)\left[i{\pi }_{lm}\left(\theta \right){\hat{\mathbf{e}}}_{\theta }+{\tau }_{lm}\left(\theta \right){\hat{\mathbf{e}}}_{\varphi }\right],$$

$${\mathbf{E}}_{\mathrm{rad}}\left(r,\theta ,\varphi \right)=\sum _{l=1}^{\infty }\sum _{m=-l}^l{E}_{lm}\left[a_{lm}{\breve{\mathbf{N}}}_{lm}+b_{lm}{\breve{\mathbf{M}}}_{lm}\right],$$

$$\begin{array}{r}{\breve{\mathbf{N}}}_{lm}=\left[{\tau }_{lm}\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }+i{\pi }_{lm}\left(\cos \theta \right){\hat{\mathbf{e}}}_{\varphi }\right]\frac{{\left[kr{h_l}^{\left(1\right)}\left(kr\right)\right]}^{\prime }}{kr}\\ +\exp \left(im\varphi \right){\hat{\mathbf{e}}}_{r}l\left(l+1\right)P_l^m\left(\cos \theta \right)\frac{{h_l}^{\left(1\right)}\left(kr\right)}{kr}\exp \left(im\varphi \right);\end{array}$$

$${\breve{\mathbf{M}}}_{lm}=\left[i{\pi }_{lm}\left(\cos \theta \right){\hat{\mathbf{e}}}_{\theta }-{\tau }_{lm}\left(\cos \theta \right){\hat{\mathbf{e}}}_{\varphi }\right]{h_l}^{\left(1\right)}\left(kr\right)\exp \left(im\varphi \right),$$

$$Q=\frac{{\tilde{\omega }}_1}{2{\tilde{\omega }}_2}.$$

Full text

Electromagnetic multipoles characterized by vector spherical harmonics constitute a complete basis for light field expansions, and are playing indispensable roles in diverse branches of optics and photonics jackson1962classical ; Bohren1983_book ; RAAB_2004__Multipole . Those multipoles have been extensively studied and their general features, such as radiation patterns, polarization distributions, far-field parities and so on, have been comprehensively exploited for various applications. A rather outstanding example of this is the burgeoning field of meta-optics largely built on electric and magnetic multipoles of different orders, the interferences among which can render enormous extra freedom for manipulating light-matter interactions in both linear and nonlinear regimes jahani_alldielectric_2016 ; KUZNETSOV_Science_optically_2016 ; SMIRNOVA_Optica_multipolar_2016 ; LIU_2018_Opt.Express_Generalized .

In spite of the great achievements relying on electromagnetic multipoles, an unfortunate situation is that rare attention has been paid to those dark directions along which there are no radiations. Those directions are naturally considered to be trivial; singularities and topological features of multipoles have neither attracted much interest nor properly investigated. This is rather puzzling, as from a mathematical viewpoint, in the far field multipolar radiations are transverse and perfectly makes an elementary case of tangent vector fields on a momentum sphere. Depending on the order of multipoles, there are a finite number (at least one) of singularities (zeros or dislocations) of the vector fields. They are isolated directions where there are no radiations and for which Poincaré indexes (winding numbers or topological invariants) can be assigned. Moreover the Poincaré-Hopf theorem, or more precisely the hairy ball theorem for this specific case, can be simply applied MILNOR_1997__Topology ; RICHESON_2012__Euler ; GBUR_2016__Singular ; BERRY_2017__HalfCentury . It is even more stunning that such a situation has been stagnant, withstanding the recent overwhelming trend of incorporating topological concepts into almost every branch of physics Hasan2010_RMP ; Qi2011_RMP ; YANG_Phys.Rev.Lett._topological_2015-1 ; HUBER_2016_Nat.Phys._Topological ; Lu2014_topological ; OZAWA_2018_ArXiv180204173 . Considering that electromagnetic multipoles broadly serve as a basic language for descriptions and explanations of various optical effects, to clarify their topological features is of apparently great significance, especially for the expanding field of topological photonics Lu2014_topological ; OZAWA_2018_ArXiv180204173 .

In this work we reapproach electromagnetic multipoles from a topological perspective, focusing on distributions of singularities and their indexes of the corresponding tangent vector fields on a $S^{2}$ momentum sphere. A comprehensive map is given, pinpointing all singularities with their indexes for multipoles of different orders. It is shown explicitly that of all multipoles the index sum over the entire momentum sphere is always $2$ (Euler characteristic of a sphere), consistent with the Poincaré-Hopf theorem MILNOR_1997__Topology ; RICHESON_2012__Euler . Based on this revelation, we reinterpret the recently observed bound states in the continuum (BICs) HSU_Nat.Rev.Mater._bound_2016 from a multipolar perspective, and trace the mechanism to overlapping radiation singularities with open radiation channels. We further discover that topological charges of the BICs demonstrated have a multipolar underlying structure, which essentially are equivalent to singularity indexes of corresponding multipolar radiations. The incorporation of topological terminology into classical multipoles can potentially rejuvenate many multipole related studies, paving the way for introducing established topological concepts and advanced functionalities into more disciplines both within and beyond photonics GBUR_2016__Singular ; BERRY_2017__HalfCentury .

Electromagnetic multipoles can be categorized into magnetic and electric groups, the radiated electric (magnetic) fields of which are characterized by vector spherical harmonics $\textbf{M}_{lm}$ ($\textbf{N}_{lm}$) and $\textbf{N}_{lm}$ ($\textbf{M}_{lm}$), respectively jackson1962classical ; Bohren1983_book ; Supplemental_Material . For analysis convenience one specific expression of $\textbf{M}_{lm}$ is chosen: Bohren1983_book ; Supplemental_Material :

[TABLE]

Here we have dropped the radially dependent terms to focus on the tangent fields; $P_{l}^{m}\left({\cos\theta}\right)$ ($-l\leq m\leq l$) denotes associated Legendre polynomials BRONSHTEIN_2007__Handbook ; basis vectors $\mathbf{\hat{e}}_{\theta}$, $\mathbf{\hat{e}}_{\phi}$ and $\mathbf{\hat{e}}_{r}$, polar and azimuthal angles $\theta$ and $\phi$ are shown in Fig. 1(a); $\textbf{M}_{lm}\cdot\mathbf{\hat{e}_{r}}=0$ and $\textbf{N}_{lm}\cdot\textbf{M}_{lm}=0$. As a first step, we investigate the poles with $|\cos\theta|=1$, where $\phi$ is not defined. Except for $|m|=1$, both are singularities and the index is Supplemental_Material :

[TABLE]

When $m=\pm 1$, the poles are not singular, for which we can still assign the index of [math] and thus Eq. (2) is still applicable Supplemental_Material .

Next we turn to regions excluding the poles, where for $m=0$ there are no singularities Supplemental_Material . When $m\neq 0$, according to Eq. (1), a singularity requires that both components along $\mathbf{\hat{e}}_{\theta}$ and $\mathbf{\hat{e}}_{\phi}$ are zero, and thus indexes and positions of singularities are Supplemental_Material :

[TABLE]

For $0<\theta<\pi$ and $0\leq\phi\leq 2\pi$: both $\sin(m\phi)=0$ and $\cos(m\phi)=0$ have $2|m|$ solutions; while equations of $P_{l}^{m}\left({\cos\theta}\right)=0$ and $\frac{dP_{l}^{m}({\cos\theta})}{d\theta}=0$ have $l-|m|$ and $l-|m|+1$ solutions, respectively BRONSHTEIN_2007__Handbook . It means that except for the poles, there are $2|m|(l-|m|)$ singularities with index of $-1$ and $2|m|(l-|m|+1)$ singularities with index of $+1$ (this is still a valid statement for $m=0$ Supplemental_Material ). Basically the index sum for all singularities is $2(1-|m|)+2|m|(l-|m|+1)-2|m|(l-|m|)$, which would always make up to $2$ and agrees with the Poincaré-Hopf theorem MILNOR_1997__Topology ; RICHESON_2012__Euler .

Up to now, we have discussed only magnetic multipoles characterized by ${{\bf{M}}_{lm}}$ with one specific expression shown in Eq. (1). For the other expression of ${{\bf{M}}_{lm}}$ [with $\sin(m\phi)$ and $\cos(m\phi)$ interchanged and $\mathbf{\hat{e}}_{\theta}$ changed to $-\mathbf{\hat{e}}_{\theta}$ Bohren1983_book ; Supplemental_Material ), identical distribution map for singularities and indexes would be obtained except for an interchange between $\sin(m\phi)$ and $\cos(m\phi)$ in Eq. (3). For electric multipoles characterized by ${{\bf{N}}_{lm}}$, the duality of Maxwell equations requires that the corresponding tangent fields can be obtained directly from those of ${{\bf{M}}_{lm}}$ (the terms with radial dependence are dropped for studies of tangent fields) by the following transformation Bohren1983_book ; Supplemental_Material : $\mathbf{\hat{e}}_{\theta}\rightarrow\mathbf{\hat{e}}_{\phi}$ and $\mathbf{\hat{e}}_{\phi}\rightarrow-\mathbf{\hat{e}}_{\theta}$. Since the transformation above would neither change the position nor the index of the singularity Supplemental_Material , as a result both electric and magnetic multipoles of the same order $(l,m)$ share the same singularity and index distributions, as summarized in Table  1. For example, Fig. 1(b) shows radiation patterns and pinpoint some representative singularities with indexes for three multipoles (see Ref. Supplemental_Material for more scenarios). The vector field patterns close to singularities of different indexes are shown in Fig. 1(c). For the case of index $+1$, the singular point could be a source (sink) or a centre, depending on which field (electric or magnetic) we employ to show the field lines Supplemental_Material .

Table 1 tells that a singularity with index larger than $+1$ is not accessible with an individual multipole. In contrast, if a series of multipoles of different orders and/or natures are combined together, a singularity with larger positive index is always accessible, which can be simply justified by the following arguments: for a specific index, we can always predesign a continuous tangent vector field that includes such a singularity; the electromagnetic multipoles constitute an orthogonal and complete bases for vector field expansion; the predesigned vector field can always be expanded into a set of multipoles, with expansion coefficients of $a_{lm}$ and $b_{lm}$ that are associated with $\textbf{N}_{lm}$ and $\textbf{M}_{lm}$, respectively Supplemental_Material .

According to the Poincaré-Hopf theorem MILNOR_1997__Topology ; RICHESON_2012__Euler , the index sum for all singularities over the momentum sphere should always be $2$, being the radiations from an isolated multipole or a set of multipoles. The simplest allowed case of this theorem is that there is only one singularity (at least one singularity) and the index thus has to be $+2$. This corresponds to (generalized) Kerker conditions (Kerker multipoles) or (generalized) Huygens sources jahani_alldielectric_2016 ; KUZNETSOV_Science_optically_2016 ; LIU_2018_Opt.Express_Generalized . The field pattern close to a singularity of index $+2$ (also termed as dipole singularity MILNOR_1997__Topology ; RICHESON_2012__Euler ) is shown in Fig. 2(a). This type of singularity can be obtained from interferences of a pair of electric and magnetic multipoles of the same order, magnitude and phase ($a_{lm}=b_{lm}$, $m=\pm 1$), but opposite parities  Liu2014_ultradirectional ; LIU_Phys.Rev.Lett._generalized_2017 , with the first three cases of overlapping dipoles, quadrupoles, and octupoles shown in Figs. 2(b)-(d).

For a tangent and continuous vector field $\mathbb{M}_{lm}$ obtained from Eq. (1) by forcing the following transformation (forbidden for $m=\pm 1$ that breaks the continuity of the vector field Supplemental_Material ) $\mathbf{\hat{e}}_{\theta}\rightarrow-\mathbf{\hat{e}}_{\theta}$ (or equivalently $\mathbf{\hat{e}}_{\phi}\rightarrow-\mathbf{\hat{e}}_{\phi}$):

[TABLE]

The singularity distribution table for $\mathbb{M}_{lm}(|m|\neq 1)$ would be the same as Table 1, except that the index in the second column is changed from $1-|m|,~{}-1,~{}+1$ to $1+|m|,~{}+1,~{}-1$ Supplemental_Material . This means that for such a transformed vector field, the index for both singularities at poles would be $1+|m|$. For example, in Fig. 2(e) we show the field patterns close to a singularity of index $+3$, which is present on the poles of the radiation pattern of $\mathbb{M}_{22}$ (shown in Fig. 2(f); despite the index difference, the radiation pattern is the same as that of $\textbf{M}_{22}$ Supplemental_Material ). As has been already explained, $\mathbb{M}_{22}$ cannot be represented by an individual multipole, but by a series of multipoles, with the multipolar composition (normalized magnitudes of expansion coefficients) shown in Fig. 2(g).

After mapping out the singularities and their indexes for electromagnetic multipoles, now we proceed to show how those properties could be powerfully useful, taking recently demonstrated BICs as an example HSU_Nat.Rev.Mater._bound_2016 . For convenience of multipolar analysis, the periodic structure can be approached from an alternative reductionist perspective LEPETIT_2014_Phys.Rev.B_Controlling ; MONTICONE_2017_NewJ.Phys._Bound ; LIU_2017_ACSPhotonics_Beam ; DOELEMAN_2018_Nat.Photonics_Experimentala ; HA_2018_Nat.Nanotechnol._Directionala : it can be treated as an infinite ensemble of radiating items (unit cells); the overall optical properties of the periodic structure can be interpreted as interferences of radiations from all the unit cells. For an isolated radiating item, there are an infinite number of open out-coupling radiation channels, corresponding to all the points (directions) on the momentum sphere. While for the periodic structure, there are only a finite number of such channels, corresponding to diffractions channels of different orders along certain directions. From this perspective, the formation of a BIC can be attributed to overlapping of the radiation singularities of each unit cell with the open diffraction channels of the periodic structure. An index could be assigned to the singularity of the unit-cell radiation that overlaps with the open channel (there can be extra singularities not overlapping with such channels), which is exactly the topological charge of the induced BIC HSU_Nat.Rev.Mater._bound_2016 ; ZHEN_2014_Phys.Rev.Lett._Topological .

To further crystallize what has been stated above, we turn to the photonic crystal slabs of square or hexagonal lattices, which are wieldy employed for investigations into BICs HSU_Nat.Rev.Mater._bound_2016 ; YANG_2014_Phys.Rev.Lett._Analytical ; ZHEN_2014_Phys.Rev.Lett._Topological ; HSU_Nature_observation_2013-1 ; GUO_Phys.Rev.Lett._topologically_2017 ; KODIGALA_Nature_lasing_2017 ; ZHANG_2018_Phys.Rev.Lett._Observation . For a direct comparison between our results from multipolar analysis with those established ones, and to connect the singularity indexes of multipolar radiations to the confirmed topological charges of BICs, we investigate structures similar to those studied in Refs. ZHEN_2014_Phys.Rev.Lett._Topological ; HSU_Nature_observation_2013-1 . The photonic crystal slabs of square or hexagonal lattices of circular air holes (slab refractive index $n=1.5$; width $w$; lattice constant $p$; air hole diameter $p/2$) are studied, which are shown as insets in Figs. 3(a) and (d). For both cases, we show the dispersion curves of two representative bands (one TE-like and one TM-like bands; $\widetilde{\omega}$ is the complex mode eigenfrequency Supplemental_Material ).

Firstly we check the symmetry-protected BICs located on the $\Gamma$-points. For both lattices, there are two such BICs on both the TE-like and TM-like bands [four BIC points A-D are indicated in Figs. 3(a) and (d), with $w=p/2$]. For each of them, we show in Figs. 3(b)-(c) and (e)-(f) the corresponding major dominant multipolar components (see Ref. Supplemental_Material for details of multipolar expansions based on near-field currents) and the far-field radiation patters. The four symmetry-protected BICs have a common feature that unit-cell radiation corresponds to that of a dominant individual multipole with fixed order and the radiation singularity overlaps with the allowed radiation channel (radiation is zero along $\mathbf{z}$ direction). This agrees with our previous discussions with regard to the multipolar interpretation of the underlying mechanism of BICs. Moreover, based on the dominant multipolar components, for each singularity along $\mathbf{z}$ direction an index can be directly assigned according to the $|\cos(\theta)|=1$ category in Table 1 (other minor multipoles not shown modifies a bit the overall radiation patterns, while not affecting singularity and index distributions Supplemental_Material ), which is exactly the topological charge revealed in Ref. ZHEN_2014_Phys.Rev.Lett._Topological .

Now we turn to the off-$\Gamma$-point accidental BICs ZHEN_2014_Phys.Rev.Lett._Topological ; HSU_Nature_observation_2013-1 , taking the indicated point E ($k_{x}p/2\pi=0.299$) in Fig. 3(a) on the TM-like band for example [see similar analysis in Ref. Supplemental_Material for indicated BIC point F in Fig. 3(d) of the hexagonal lattice]. We show respectively both its multipolar composition and reconstructed far-field patterns of the unit cell in Figs. 4(d) and (b). As expected, Fig. 4(b) shows that a pair of singularities (both of index $+1$ guaranteed by $\sigma_{z}$ mirror symmetry) overlap with the allowed diffraction directions [$\phi_{d}=0$ and $\theta_{d}=\arcsin(k_{x}/k_{0})$, where $k_{0}$ is total the angular wavenumber]. Moreover, the index $+1$ of the singularity agrees with the topological charge previously revealed ZHEN_2014_Phys.Rev.Lett._Topological . This BIC would move to another position E′ ($k_{x}p/2\pi=0.377$) by changing $w$, as is shown in Fig. 4(a). It is clear from Figs. 4(c) and (d) that the multipolar composition changes accordingly such that the multipolar singularity pair coincides with the new diffraction direction, which ensures that the BIC can move smoothly to the new position of different $k_{x}$ (see Ref. Supplemental_Material for decay states which correspond to separated singularities and open channels).

Also according to Fig. 4(a), the symmetry-protected BIC is fixed at $\Gamma$-point with changing $w$, which can certainly be explained through symmetry analysis ZHEN_2014_Phys.Rev.Lett._Topological ; KODIGALA_Nature_lasing_2017 ; ZHANG_2018_Phys.Rev.Lett._Observation ; SAKODA_2004__Optical . Our multipolar interpretation gives a different insight: for this type of BIC, the unit-cell radiation is represented by multipoles of order $|m|\neq 1$ Supplemental_Material , the radiations of which are intrinsically zero along $\mathbf{z}$ direction (see Figs. 1 and  3, and Ref. Supplemental_Material ) and can thus maintain the BIC on $\Gamma$-point. In sharp contrast, for the accidental BICs shown in Fig. 4, zero radiation is accidental that originates from completely destructive multipolar interferences. Changing geometric parameters would alter the multipolar ratios and thus also the singularity positions, inevitably shifting the BICs. Moreover, as we have argued before, if the unit-cell radiation is represented by a dominant individual multipole, it is not possible to obtain positive index (topological charge) larger than $+1$ (see Table  1), no matter which band we rely on. Though symmetry analysis reveals that large positive index is not forbidden  ZHEN_2014_Phys.Rev.Lett._Topological ; KODIGALA_Nature_lasing_2017 ; ZHANG_2018_Phys.Rev.Lett._Observation ; SAKODA_2004__Optical , here we make a further step to show more microscopically that, to obtain a larger positive index, at least two combined multipoles are required (see Fig. 2). This is probably synonymous with the phenomenon that: in previous studies large negative-index BIC have been widely achieved ZHEN_2014_Phys.Rev.Lett._Topological ; KODIGALA_Nature_lasing_2017 ; ZHANG_2018_Phys.Rev.Lett._Observation ; SAKODA_2004__Optical , whereas their positive counterparts rarely manifest themselves.

In conclusion, we revisit electromagnetic multipoles from a topological perspective, and provide a comprehensive map for the distributions of singularities and their indexes for multipolar radiations. It’s shown that for an individual multipole, it is not possible to obtain singularity index larger than +1; while for the combination of a series of multipoles there is no limit for its singularity, but only requires that index sum of all singularities over the momentum space has to be $2$. Singularities of multipolar radiations are synonymous with the formation of BICs, as long as they coincide with open radiation channels. Based on this multipolar revelation, we further uncover the subtle equivalence between singularity indexes of multipoles and topological charges of BICs. Considering the ubiquitous roles of electromagnetic multipoles all across photonics, our work will accelerate the pervasion of topological concepts into more optical branches, bring unperceived opportunities for various applications. Furthermore, multipoles basically serve as a fundamental tool and language for many other fields involving wave effects, on which our work generally shed new light from a fundamental topological perspective.

We acknowledge the financial support from National Natural Science Foundation of China (Grant No. 11874026, 11404403 and 11874426), and the Outstanding Young Researcher Scheme of National University of Defense Technology.

I. Expressions for vector spherical harmonics

Vector spherical harmonics can be generated from a scalar function $\Psi$ through the following relations Bohren1983_book :

[TABLE]

where $\mathbf{r}$ is the radial vector and $k$ is angular wavenumber. $\Psi$ can be an even or odd function with respect to azimuthal angle $\phi$, which are respectively expressed as Bohren1983_book :

[TABLE]

where $P_{l}^{m}\left({\cos\theta}\right)$ ($-l\leq m\leq l$) are associated Legendre polynomials and $z_{n}(kr)$ is the spherical Bessel or Hankel function BRONSHTEIN_2007__Handbook . Eqs. (S1) and (S2) would lead to the two sets of expressions for both $\mathbf{M}$ and $\mathbf{N}$ Bohren1983_book (in this work as we fucus on the radiated tangent vector fields, the $\mathbf{r}$-dependent terms are dropped and all components along $\mathbf{r}$ are neglected):

[TABLE]

It is clear from Eq. (S3) that we can convert ${{\bf{M}}^{e,o}_{lm}}$ directly to ${{\bf{N}}^{e,o}_{lm}}$ through the following transformation (duality transformation): $\mathbf{\hat{e}}_{\theta}\rightarrow\mathbf{\hat{e}}_{\phi}$ and $\mathbf{\hat{e}}_{\phi}\rightarrow-\mathbf{\hat{e}}_{\theta}$, which is guaranteed by the duality of Maxwell equations jackson1962classical . In the main letter, we have chosen ${{\bf{M}}^{e}_{lm}}$ for analysis. The analysis shown below about singularity and index distributions can be easily extended to ${{\bf{M}}^{o}_{lm}}$ and ${{\bf{N}}^{e,o}_{lm}}$, and the same conclusions can be drawn.

II. Positions and indexes for singularities of the multipolar radiations

For the multipole characterized by ${{\bf{M}}^{e}_{lm}}$, the singularity requires that

[TABLE]

The Poincaré index of the singularity located at $(\theta_{s},\phi_{s})$ can be defined as RUDIN :

[TABLE]

where $\mathbf{\gamma}$ is a continuous differentiable closed contour containing this singularity (all other singularities locate outside of the contour) and $d\alpha$ is the rotation angle of the vector field with respect to ${\mathbf{\hat{e}}_{r}}(\theta_{s},\phi_{s})$. Basically when the index is positive (negative), the rotation of the field vector is in the same (opposite) direction that the contour is traversed.

II.I Singularity and index at poles of $\theta=0$ or $\pi$

$\blacksquare$   For the case of $m\neq\pm 1$

It is easy to see that $|{{\bf{M}}^{e}_{lm}}|=0$ when $|m|\neq 1$ and $\theta=0$ or $\pi$ BRONSHTEIN_2007__Handbook , which means that multipoles with $|m|\neq 1$ would have two singularities located on the two poles. For $\theta=0$ or $\pi$, the total rotation angle of the vector field close to singularity (along $\mathbf{z}$ or $-\mathbf{z}$) can be expressed as (here we assume that the contour is traversed in a counter-clockwise manner):

[TABLE]

where $\sigma$ is the helicity of the vector field in the rotating orthonormal frame of $({\mathbf{\hat{e}}_{\theta}},~{}{\mathbf{\hat{e}}_{\phi}},~{}{\mathbf{\hat{e}}_{r}})$. Similar to the definition of circularly polarized light Yariv2006_book_photonics : $\sigma=1$ when the phase lag between the ${\mathbf{\hat{e}}_{\theta}}$ and ${\mathbf{\hat{e}}_{\phi}}$ components is $\beta_{{\mathbf{\hat{e}}_{\phi}}}-\beta_{{\mathbf{\hat{e}}_{\theta}}}=-\pi/2$ and thus the vector field is rotating in a counter-clockwise manner in the rotating frame; or $\sigma=-1$ when $\beta_{{\mathbf{\hat{e}}_{\phi}}}-\beta_{{\mathbf{\hat{e}}_{\theta}}}=\pi/2$ and the vector field is rotating in a clockwise manner in the rotating frame; an extra special case is that $m=0$, which means that there is only ${{\mathbf{\hat{e}}_{\phi}}}$ or ${{\mathbf{\hat{e}}_{\theta}}}$ component [see Eq. (S3)], indicating a linear polarization and thus $\sigma=0$. Now the origin of the two terms on the right-hand side of Eq. (S6) is clear: the term $2\pi$ comes from the rotation of the frame itself; while the term $2\pi|m|\sigma$ comes from the ration of the vector field within the rotating frame.

For all the vector fields expressed by Eq. (S3), it is easy to confirm that $\sigma=0$ for $m=0$ and $\sigma=-1$ for $|m|\neq 0$. From Eqs. (S5) and (S6) we obtain:

[TABLE]

which is valid for $m=0$. Equation (S7) confirms that for multipoles characterized by both ${{\bf{M}}_{lm}}$ and ${{\bf{N}}_{lm}}$, the singularity index at the poles are identical, confirming that the duality transformation would neither change the position nor the index of the singularity.

$\blacksquare$   For the case of $m=\pm 1$

When $m=\pm 1$, $|{{\mathbf{M}}^{e}_{lm}}|\neq 0$ for $\theta=0$ or $\pi$. This means that the vector field is not singular at the south or north poles, for which we can assign an index of [math], and thus Eq. (S7) is still valid.

$\blacksquare$   Singularity and index distribution for vector fields of $\mathbb{M}^{e}_{lm}$ at $\theta=0$ or $\pi$.

For a tangent and continuous vector field $\mathbb{M}^{e}_{lm}$ obtained from $\mathbf{M}^{e}_{lm}$ by forcing the following transformation $\mathbf{\hat{e}}_{\theta}\rightarrow-\mathbf{\hat{e}}_{\theta}$ (or equivalently $\mathbf{\hat{e}}_{\phi}\rightarrow-\mathbf{\hat{e}}_{\phi}$):

[TABLE]

According to the previous analysis, the transformation above has flipped the helicity from $\sigma$ to $-\sigma$ within the rotating frame. As a result, the index for $\mathbb{M}^{e}_{lm}$ can be expressed as:

[TABLE]

As will be explained below, Eq. (S9) is only valid for $m\neq\pm 1$.

$\blacksquare$   Why a continuous $\mathbb{M}^{e}_{lm}$ Does not exist for $m=\pm 1$.

When $m=\pm 1$, similar to $\mathbf{M}^{e}_{lm}$, $|{{\mathbb{M}}^{e}_{lm}}|\neq 0$ for $\theta=0$ or $\pi$. This means $\mathbb{M}^{e}_{l,\pm 1}$ would be non-vanishing at the south or north poles. For a non-vanishing vector field at the poles, the continuity requires that all field vectors point to the same direction close to the poles. For $\mathbf{M}^{e}_{l,\pm 1}$, this is exactly the case. Unfortunately this is not the case anymore for $\mathbb{M}^{e}_{l,\pm 1}$, as according to our previous analysis, the field vectors close to pole would be exactly the same as that of a singularity of index $1+|m|$, rotating by an angle of $2(1+|m|)$ along the closed contour. For example, in Fig. S1 we show the vector field patterns around the poles for both $\mathbf{M}^{e}_{l1}$ and $\mathbb{M}^{e}_{11}$. It is clear that $\mathbf{M}^{e}_{11}$ can be both non-vanishing and continuous across the poles while this is not possible for $\mathbb{M}^{e}_{11}$.

II.II Singularity and index when $\theta\neq 0$ or $\pi$

In the regions excluding the poles, azimuthal angle $\phi$ is well defined and a singularity requires that both the both components along ${{\mathbf{\hat{e}}_{\phi}}}$ and ${{\mathbf{\hat{e}}_{\theta}}}$ vanish.

$\blacksquare$   Singularity and index for $\mathbf{M}^{e}_{lm}$.

For $\mathbf{M}^{e}_{lm}$, a singularity requires that:

[TABLE]

or

[TABLE]

For both scenarios shown in Eqs. (S10) and  (S11), the common feature is that the singularity comes from the crossings of two nodal lines, which means that the index of the singularity would be either $+1$ (sink, source or center singularity) or $-1$ (saddle singularity). The exception of this is that when $m=0$: Eq. (S10) does not have solution; while for Eq. (S11), the nodal line of $\sin(m\phi)=0$ would then evolve into a nodal plane, for which we can assign an index of [math]. It is worth mentioning that for our study into multipolar radiations, center-type singularity are equivalent to sink (source)-type singularity, since if it is a sink (source) for electric field, then it becomes a center for magnetic field, and vice versa.

For $0<\theta<\pi$ and $0\leq\phi\leq 2\pi$: both equations of $\sin(m\phi)=0$ and $\cos(m\phi)=0$ have $2|m|$ solutions ($m\neq 0$); while equations of $P_{l}^{m}\left({\cos\theta}\right)=0$ and $\frac{dP_{l}^{m}({\cos\theta})}{d\theta}=0$ have $l-|m|$ and $l-|m|+1$ solutions ($m\neq 0$), respectively BRONSHTEIN_2007__Handbook . As a result, Eq. (S10) corresponds to $2|m|(l-|m|)$ singularities while Eq. (S11) corresponds to $2|m|(l-|m|+1)$ singularities. Both expressions are applicable for the special case of $m=0$, when there are no singularities located between the two poles. Rotational symmetry guarantees that the singularities described by the same equation share the same index of $-1$ or $+1$. When $m\neq 0$, simple vectorial analysis of the two nodal lines (schematically shown in Fig. S2) tells us that the singularities described in Eq. (S10) are of saddle-type and show index of $-1$; while those described in Eq. (S11) are of source (sink or centre)-type and thus show index of $+1$. This can be summarized as follows ($m\neq 0$):

[TABLE]

Over the entire momentum sphere, the index sum for all singularities is:

[TABLE]

which would always make up to $2$, agreeing with the Poincaré-Hopf theorem MILNOR_1997__Topology ; RICHESON_2012__Euler .

$\blacksquare$   Singularity and index for \mathbb{M}^{e}_{lm}$$(|m|\neq 1).

During the transformation of $\mathbf{\hat{e}}_{\theta}\rightarrow-\mathbf{\hat{e}}_{\theta}$ (or equivalently $\mathbf{\hat{e}}_{\phi}\rightarrow-\mathbf{\hat{e}}_{\phi}$), the positions of the singularities would be fixed [see Eqs. (S10) and (S11)] while the sign of the index would be flipped: a saddle would be transformed to a sink (source) and vice versa. As a result, the singularity and index distributions for $\mathbb{M}^{e}_{lm}(|m|\neq 1)$ is as follows ($m\neq 0$):

[TABLE]

For $\mathbb{M}^{e}_{lm}(|m|\neq 1)$, over the entire momentum sphere, the index sum for all singularities is:

[TABLE]

which would also always make up to $2$, as is required by Poincaré-Hopf theorem MILNOR_1997__Topology ; RICHESON_2012__Euler . The singularities and their index distributions for $\mathbb{M}^{e}_{lm}$ $(|m|\neq 1)$ has been summarized in Table  S1.

III. Radiation patterns, singularities and indexes for $\mathbf{M}^{e}_{lm}$ with $l=1\rightarrow 3$ and $m=0\rightarrow l$

Figure S3 shows the radiation patterns for multipoles described by $\mathbf{M}^{e}_{lm}$ with $l=1\rightarrow 3$ and $m=0\rightarrow l$. Singularities that represent each category in Table I in the main letter are pinpointed and the corresponding indexes are specified.

IV. Multipolar expansions and the radiation patterns for unit cells of infinite periodic structures

For finite particles or particles clusters, the radiated fields along all directions can be (numerically) calculated, which can then be directly expanded to spherical harmonics of different orders with coefficients of $a_{lm}$ and $b_{lm}$. Then the magnitudes and phases of all electric and magnetic multipoles would be obtained, as has been finely implemented in most previous studies jahani_alldielectric_2016 ; KUZNETSOV_Science_optically_2016 ; LIU_2018_Opt.Express_Generalized . This approach of radiation expansion (expansion of the radiated fields) can not be directly implemented for unit cells of infinite periodic structures, as the radiation of infinite periodic structures is only known and allowed along certain directions of diffractions.

For infinite periodic structures, we can carry out the multipolar expansion based on the near-field polarization currents $\bf{J(r)}$ of the unit cell. The multipolar coefficients of $a_{lm}$ and $b_{lm}$ can be obtained through the following relations (when there are no effective magnetization currents) jackson1962classical ; GRAHN_NewJ.Phys._electromagnetic_2012 :

[TABLE]

where $\eta$ is the vacuum impendence; ${Q_{lm}}={1\over{{{\left[{l\left({l+1}\right)}\right]}^{1/2}}}}{\left[{{{2l+1}\over{4\pi}}{{\left({l-m}\right)!}\over{\left({l+m}\right)!}}}\right]^{1/2}}$; the integration is conducted within the unit-cell region V. As long as the electric field distributions ${\bf{E}}\left({\bf{r}}\right)$ within the unit cell is obtained, the corresponding electric polarzation currents can be obtained as:

[TABLE]

where ${\varepsilon_{r}}\left({\bf{r}}\right)$ denotes the permittivity distribution and ${\varepsilon_{b}}$ is the permittivity of the background medium. Vector functions $\mathbf{S}_{lm}$ and $\mathbf{T}_{lm}$ are:

[TABLE]

for which ${\tau_{lm}}\left(\theta\right)={d\over{d\theta}}P_{l}^{m}\left({\cos\theta}\right)$, ${\pi_{lm}}\left(\theta\right)={m\over{\sin\theta}}P_{l}^{m}\left({\cos\theta}\right)$, and ${\Pi_{l}}\left({kr}\right)=kr{j_{l}}\left({kr}\right)$ is Riccati-Bessel function [$j_{l}(kr)$ is the spherical Bessel function of the first kind] BRONSHTEIN_2007__Handbook .

With the expansion coefficients $a_{lm}$ and $b_{lm}$ obtained through Eqs. (S18) and (S19), the radiated fields from the unit cell can be reconstructed as follows:

[TABLE]

where ${E_{lm}}={{{i^{l+1}}\left({2l+1}\right)}\over 2}\sqrt{{{\left({l-m}\right)!}\over{l\left({l+1}\right)\left({l+m}\right)!}}}$; and vector functions $\bf{\breve{N}}$lm and $\bf{\breve{M}}$lm are:

[TABLE]

for which ${h_{l}}^{\left(1\right)}\left({kr}\right)$ is the spherical Hankel function of the first kind BRONSHTEIN_2007__Handbook .

V. Numerical calculations of the dispersion bands, the mode Q-factors and the field distributions within the unit cells of photonic crystal slabs

For the photonic crystal slabs studied in this work, we employ a commercial software package COMSOL Multiphysics(https://www.comsol.com) to calculate the the dispersion bands. Each mode located on the band is characterized by complex eigenfrequencies with both real and imaginary parts: $\widetilde{\omega}=\widetilde{\omega}_{1}+i\widetilde{\omega}_{2}$. With $\widetilde{\omega}$ obtained, the Q-factor of the mode can be calculated as:

[TABLE]

Both the complex eigenfrequencies and electric fields within the unit cell ${\bf{E}}\left({\bf{r}}\right)$ that are employed to obtain the electric polarization currents ${\bf{J}}\left({\bf{r}}\right)$ [see Eq. (S20)], can be calculated using COMSOL.

VI. Decay (non-BIC) states and separations of multipolar singularities from open radiation channels

As is shown in the main letter, the formation of decay-free BICs can be attributed to the overlapping of multipolar singularity with the open radiation channel. When the singularity is separated from the open channel and radiation is not zero along the diffraction direction, energy leakage is allowed and then the BICs would evolve into decay states. In Fig. S4(a) [the same as Fig. 4(a) in the main letter] we have indicated two decay states by G ($\phi_{d}=0$ and $\theta_{d}\approx 16.6^{\circ}$) and H ($\phi_{d}=0$ and $\theta_{d}\approx 18.7^{\circ}$), and show the corresponding radiation patterns, singularities, indexes, diffraction directions and multipolar compositions in Figs. S4(b)-(d). For both cases, the singularities (G: $\phi_{s}=0$ and $\theta_{s}\approx 23.5^{\circ}$; H: $\phi_{s}=0$ and $\theta_{s}\approx 26.5^{\circ}$) do not overlap with the open radiation channels to prevent energy leakage, and thus the indicated modes are not bound states, showing finite Q-factors.

VII. Multipolar compositions for the $\Gamma$-point BICs

Specific multipolar compositions (normalized magnitudes of expansion coefficients) for the $\Gamma$-point BICs indicated in the main letter by A-D and $\mathbf{B^{\prime}}$ are summarized in Fig. S5. For all cases, the multipoles of $|m|=1$ (the radiations of which are not zero at poles) do not exist and thus the radiations are intrinsically zero and singular along $\mathbf{z}$ direction.

VIII. Off-$\Gamma$-point BICs of the hexagonal photonic crystal slab

Besides the square photonic crystal slab, off-$\Gamma$-point accidental BICs are also supported by the hexagonal photonic crystal slab, one of which is indicated by F ($k_{\|}p/2\pi=0.364$) in Fig. S6(a) [geometric parameters are the same as those in Fig. 3(d) in the main letter and this indicated point is on the TM-like branch along $\Gamma$-M direction]. This bound state would moved to another position F′ ($k_{\|}p/2\pi=0.438$) by changing $w$ from $p/2$ to $p$, as is shown in Fig. S6(a). We further show respectively the multipolar compositions and reconstructed far-field patterns in Figs. S6(b)-(d). Similar to the square photonic crystal slab, for both BICs the singularities of index +1 overlaps with the allowed diffraction directions (F: $\phi_{d}=30^{\circ}$ and $\theta_{d}\approx 31.7^{\circ}$; F′: $\phi_{d}=30^{\circ}$ and $\theta_{d}\approx 50^{\circ}$), where $\phi_{d}=\arccos(k_{x}/k_{\|})$, $\theta_{d}=\arcsin(k_{\|}/k_{0})$, and $k_{0}$ is total the angular wavenumber. This confirms our multipolar reinterpretation of the building mechanism of BICs and the connections between singularity indexes and their topological charges.