Selection Rules for Quasi-Bound States in the Continuum

TL;DR

This paper provides a comprehensive set of rules for controlling symmetry-protected bound states in the continuum within photonic crystal slabs, enabling precise optical resonance tuning with minimal device size.

Contribution

It offers an exhaustive catalog of selection rules for symmetry-protected BICs in various 2D photonic crystal lattices, facilitating practical device design.

Findings

Catalog of symmetry-breaking effects on BICs in six lattice types

Inverse squared dependence of quality factor on perturbation magnitude

Simplified fabrication approach for high-Q optical resonances

Abstract

Photonic crystal slabs (PCSs) are a well-studied class of devices known to support optical Fano resonances for light normally incident to the slab, useful for narrowband filters, modulators, and nonlinear photonic devices. In shallow-etched PCSs the linewidth of the resonances is easily controlled by tuning the etching depth. This design strength comes at the cost of large device footprint due to the poor in-plane localization of optical energy. In fully-etched PCSs realized in high index contrast material systems, the in-plane localization is greatly improved, but the command over linewidth suffers. This disadvantage in fully-etched PCSs, also known as high contrast gratings (HCGs), can be overcome by accessing symmetry-protected Bound States in the Continuum (BICs). By perturbing an HCG, the BIC may be excited from the free space with quality factor showing an inverse squared dependence on the magnitude of the perturbation, while inheriting the excellent in-plane localization of their unperturbed counterparts. Here, we report an exhaustive catalog of the selection rules (if and to which free space polarization coupling occurs) of symmetry-protected BICs controlled by in-plane symmetry breaking in six types of two-dimensional PCS lattices. The chosen lattices allow access to the three highest symmetry mode classes of unperturbed square and hexagonal PCSs. The restriction to in-plane symmetry breaking allows for manufacturing devices with simple lithographic fabrication techniques in comparison to out-of-plane symmetry breaking, useful for practical applications. The approach reported provides a high-level roadmap for designing PCSs supporting controllable sharp spectral features with minimal device footprints using a mature fabrication platform.

Figures

Equations

$$Q=C/{\delta }^2$$

$$Q=\frac{{\omega }_0}{d\omega }\le \frac{{\omega }_0}{\Delta \omega }=\frac{{\omega }_0}{b\Delta k^2}.$$

$${\psi }_{L,S}^{m,n},$$

$${\gamma }_e\propto {\iint }_A[{\mathbf{E}}_{inc}^{\ast }\times {\mathbf{H}}_{mode}+{\mathbf{E}}_{mode}\times {\mathbf{H}}_{inc}^{\ast }]\cdot \hat{\mathbf{z}}dxdy,$$

$${\mathbf{E}}_{inc}^{\ast }=\left[\begin{array}{c}{A}_x\\ A_y\\ 0\end{array}\right],$$

$${\mathbf{H}}_{inc}^{\ast }=\frac{1}{{\eta }_0}\left[\begin{array}{c}{A}_y\\ A_x\\ 0\end{array}\right],$$

$${\mathbf{E}}_{mode}=e^{i(\beta z-\omega t)}\left[\begin{array}{c}{E}_x\\ E_y\\ E_z\end{array}\right],$$

$${\mathbf{H}}_{mode}=e^{i(\beta z-\omega t)}\left[\begin{array}{c}{H}_x\\ H_y\\ H_z\end{array}\right],$$

$${\gamma }_e\propto {\iint }_A[A_x({\eta }_0{H}_y+E_x)+A_y({\eta }_0{H}_x+E_y)]dxdy.$$

$$\displaystyle\gamma_{e}\propto\iint_{A}\big{[}A_{x}(c_{1}\partial_{x}E_{z}+c_{2}\partial_{y}H_{z})+$$

$$\displaystyle A_{y}(c_{1}\partial_{y}E_{z}+c_{2}\partial_{x}H_{z})\big{]}\,dx\,dy$$

$$c_1=\frac{1}{i\beta }\frac{1+{\epsilon }_r(x,y)}{1-{\epsilon }_r(x,y)}$$

$$c_2=\frac{2{\eta }_0}{i\beta }\frac{1}{1-{\epsilon }_r(x,y)}$$

$${\gamma }_e\propto ⟨\left[\begin{array}{cc}{A}_x& A_y\end{array}\right]\cdot \left[\begin{array}{cc}{c}_1{\partial }_x& c_2{\partial }_y\\ c_1{\partial }_y& c_2{\partial }_x\end{array}\right]\cdot \left[\begin{array}{c}{E}_z\\ H_z\end{array}\right]⟩,$$

$${\gamma }_e\propto {\iint }_A{c}_j{\partial }_i\psi dxdy$$

$$H^0{\psi }^0={\mathcal{E}}^0{\psi }^0$$

$$H\psi =\mathcal{E}\psi$$

$${\psi }_n^1=\sum _{m\ne n}\frac{⟨({\psi }_m^0{)}^{\ast }V{\psi }_n^0⟩}{{\mathcal{E}}_m-{\mathcal{E}}_n}{\psi }_m^0.$$

$${\Gamma }_{{\psi }^1}={\Gamma }_V\otimes {\Gamma }_{{\psi }^0}.$$

$${\gamma }_e\propto {\iint }_A{c}_j{\partial }_i\psi dxdy={\iint }_A{c}_j{\partial }_i({\psi }^0+{\psi }^1)dxdy$$

$$c_{H^0}=\frac{2{\eta }_0}{i\beta }\frac{1}{1-{\epsilon }_r^0(x,y)},$$

$$c_V=\frac{\delta {\epsilon }_r(x,y)}{1-{\epsilon }_r^0(x,y)}{c}_{H^0},$$

$${\gamma }_e\propto {\iint }_A{c}_{H^0}{\partial }_i{\psi }^{0}dxdy+{\iint }_A{c}_{H^0}{\partial }_i{\psi }^{1}dxdy$$

$$+{\iint }_A{c}_V{\partial }_i{\psi }^{0}dxdy+{\iint }_A{c}_V{\partial }_i{\psi }^{1}dxdy.$$

$${\Gamma }_{{\partial }_i}={\Gamma }_V\otimes {\Gamma }_{{\psi }^0}.$$

$${\gamma }_V\propto {\iint }_A{\partial }_i(V{\psi }^0)dxdy$$

$${\iint }_A{c}_{H^0}{\partial }_i{\psi }^{0}dxdy,$$

$$FoM=\left|\iint {\chi }^{(3)}(x,y)E^{\ast }({\omega }_1)E({\omega }_2)E^{\ast }({\omega }_3)E^{\ast }({\omega }_4)dxdy\right|$$

$$FoM=\left|{\chi }_{Si}^{(3)}{F}_3^{\ast }{F}_4^{\ast }{\iint }_A{F}_y(x,y)dxdy\right|$$

$$\varphi =2\alpha .$$

$${\Phi }_{LCP}^r$$

Full text

I Introduction

Enhancement of light-matter interactions is a key capability for improving and expanding the functionality of a wide gamut of photonic devices. Spatially and temporally confining light enables compact planar optical modulators with fast switching speeds Xu et al. (2007); Nguyen et al. (2012); Melikyan et al. (2014); Gao et al. (2014), narrowband bandpass filters Wang et al. (1990); Wang and Magnusson (1993); Tibuleac and Magnusson (1997); Rosenblatt et al. (1997); Cui et al. (2015), sensitive biological and refractive index sensors Zhuo and Cunningham (2015); Scullion et al. (2011); Tittl et al. (2018), efficient optical micro-electromechanical devices Antoni et al. (2011); Kemiktarak et al. (2014), novel lasers Bigelow et al. (2003); Zhang et al. (2006); Noda (2010); Chase et al. (2010); Nguyen et al. (2018a), and enhanced nonlinear Corcoran et al. (2009); Galli et al. (2010); Lin et al. (2016); Lee et al. (2014) and quantum optical phenomena Yoshie et al. (2004); Wang et al. (2018). This is conventionally achieved by the introduction of an optical cavity, which circulates optical energy, affording a photon many passes through a material.

Planar diffractive optics enable uniquely compact optical confinement in lightweight quasi-two-dimensional systems fabricated by mature micro- and nano-fabrication technologies. Traditional plasmonic materials such as gold enable strong light-matter interaction in metasurfaces Yu et al. (2011); Luk’yanchuk et al. (2010); Cui et al. (2012); Khanikaev et al. (2013); Li et al. (2015); Lee et al. (2014), but are incompatible with standard complementary metal-oxide semiconductor (CMOS) foundries. Alternative plasmonic materials are an active area of studyNaik et al. (2013); Wang et al. (2017), but without exception introduce substantial optical losses that reduce the efficiency of a photonic device. These limitations motivate exploring methods of confining optical energy without metals, restricting the optical materials to common dielectric materials such as silicon and its oxide.

A classic example of a dielectric diffractive optical element with enhanced light-matter interactions is the low-contrast grating (LCG), or guided mode resonance filter Wang et al. (1990); Wang and Magnusson (1993); Tibuleac and Magnusson (1997); Rosenblatt et al. (1997). By periodically corrugating a thin slab with subwavelength periodicity, a laterally propagating waveguide mode supported by the slab may couple to normally incident light. The leakage out of the slab interferes with the direct optical pathways (here, the Fabry-Perot resonance), producing a well-known Fano resonance Fano (1935, 1961); Limonov et al. (2017). Related phenomena have been studied for over a century, beginning with Wood’s anomalies Wood (1902); Hessel and Oliner (1965); Maystre (2012). In an LCG, the degree of corrugation can be easily controlled experimentally, and is a design parameter that directly controls the linewidth of the resonant spectral feature. In particular, for small corrugation the Q-factor of the resonance is known to be inversely related to the depth of the corrugation Fan et al. (2003). However, this attractive design feature comes with an inherent drawback: the long optical lifetime comes from the long distance the guided mode travels within the device before coupling back to free space; the device therefore needs to be of a lateral size comparable to this characteristic travel distance in order to observe a narrow spectral feature. In other words, LCGs are constrained by a tradeoff between spatial confinement (device size) and temporal confinement (Q-factor).

Another well-studied diffractive optical element is the high contrast grating (HCG) Chang-Hasnain and Yang (2012); Karagodsky et al. (2010); Chang-Hasnain (2011), known to enable compact devices due to large in-plane Bragg reflection laterally confining optical energy. Since the corrugation is deep (and, typically, complete) in HCGs, the ease of control of the Q-factor by the method present in LCGs is lost. HCGs are best-known for their broadband spectral features for this reason. However, HCGs are also known to support sharp spectral features in the form of Fano resonances Chang-Hasnain (2011); Yoon et al. (2016); Wang et al. (2016). In particular, for certain combinations of optical materials, geometries, wavelength, angle, and polarization, the Q-factor may become infinite, a phenomenon known as a “bound state in the continuum” (BIC) Plotnik et al. (2011); Hsu et al. (2013); Bulgakov and Sadreev (2014); Hsu et al. (2016). Operating near a BIC in the relevant multi-dimensional parameter space allows tuning of a resonance with finite Q-factor. Unfortunately, because of the complex and sensitive dependence on many parameters simultaneously, this control is not robust in comparison to the control in an LCG.

However, HCGs can support two classes of BICs: those excluded from coupling to free space due to symmetry constraints (or “symmetry-protected”), and those excluded for reasons unrelated to symmetry (or “accidental” Sadrieva et al. (2019)). It has been argued recently Overvig et al. (2018) that symmetry-protected BICs in HCGs are better suited than accidental BICs for creating compact optical devices with sharp spectral features. It is well-known that by reducing the symmetry Lemarchand et al. (1999); Bingham et al. (2008); Zeng et al. (2015); Nguyen et al. (2017); Liang et al. (2014); Cui et al. (2016); Qiu et al. (2012); Cui et al. (2012); Kilic et al. (2008); Foley et al. (2014); Fedotov et al. (2007); Nguyen et al. (2018b); Cueff et al. (2018) of an HCG or PCS, symmetry-protected BICs become quasi-bound in the continuum, at which point they are referred to as “quasi-BICs”Koshelev et al. (2018). Quasi-BICs couple to light at normal incidence with optical lifetimes controlled by the magnitude of the perturbation that breaks the symmetry protecting them, thereby restoring a robust design paradigm for controlling the Q-factor of a sharp spectral feature. Furthermore, it has also recently been shown Nguyen et al. (2018b); Cueff et al. (2018) that proper perturbation (including breaking vertical symmetry) allows excellent control of the band structure. Therefore, a symmetry-broken HCG inherits the benefits of both LCGs and HCGs relevant to sharp spectral features in compact devices.

In particular, a period doubling perturbation (a dimerization of an HCG) allows modes previously bound (under the light line at the edge of the first Brillouin zone (FBZ)) to be brought into the continuum, coupling to a range of angles near normal incidence to a degree controlled solely by the perturbation Wu et al. (2014); Liang et al. (2014); Zeng et al. (2015); Lan et al. (2017); Nguyen et al. (2017); Qiu et al. (2012); Nguyen et al. (2018b); Cueff et al. (2018). Consequently, a “dimerized high contrast grating” (DHCG Overvig et al. (2018)) is an excellent candidate platform for planar optical devices with both spatial and temporal confinement of light. Much of the study of DHCGs has focused on simple, one-dimensional devices, enabling control of the mode in one in-plane direction, but not in the orthogonal direction. Two-dimensional, high-index contrast PCSs with periodic perturbations are the natural extension of DHCGs that solve this limitation, and are the subject of this paper. The number of symmetries in a two-dimensional PCS is significantly greater than the simple one-dimensional case; the wealth of modal interactions between free space and two-dimensional PCSs with periodic perturbation therefore requires detailed exploration.

Symmetry-protected BICs are commonly studied in monatomic PCSs, where even/odd symmetry conditions may preclude coupling to free space at normal incidence. The BICs in diatomic PCSs (e.g., DHCGs) are subject to the analogous even/odd symmetry conditions, so that once the period doubling has folded the bound modes into the continuum, they may still be left bound in the continuum. If the relevant even/odd symmetry is broken, the symmetry-protected BICs become quasi-BICs. Both the monatomic and diatomic approaches fall under the same category of symmetry-protected BICs, but access distinct high symmetry modes (that is, modes with distinct periodicity and field profiles). Therefore, to fully utilize the available perturbations and modes, we study both monatomic and multi-atomic PCSs.

We note that the simplest method of breaking the relevant in-plane symmetries is to excite quasi-BICs with light at an incident angle just off the substrate’s normal. However, as argued in  Overvig et al. (2018), this approach is restrictive in comparison to breaking the symmetry by perturbing the structure. In particular, this method works only for light with a limited angular spread of optical momentum centered at a precise angle (to get the target Q-factor), therefore requiring devices with large footprints and precise tilt of the device relative to the source optics.

Lastly, while in this Article we describe quasi-BICs as supported by all-dielectric PCSs, we note that the Group Theory approach taken here is valid irrespective of materials system, so long as the materials are isotropic. For instance, arrays of Silicon pillars are treated the same as arrays of holes in a Silicon (or even metal) slab. The choice of materials system may be made based on the needs of the application; an analogous set of modes (as classified by their in-plane symmetries) obeying the selection rules derived here will exist, albeit with widely varying field profiles and resonant frequencies.

In this paper, we study the optical response due to in-plane perturbations applied to high-symmetry PCS lattices. Throughout, we assume these PCSs have cross-sections invariant in the out-of-plane direction, and we assume any out-of-plane symmetry introduced by the presence substrate is negligible. In Sec. II.1, we review quasi-BICs and how they spatially and temporally confine light. In Section II.2, we review the classification in the language of Group Theory of the three types of high symmetry modes supported by each of square and hexagonal PCSs (six types of modes in total). Section III.1 explores six classes of perturbed lattices chosen to target the six high symmetry modes. To determine the impact of these perturbations, Sec. III.2 derives the symmetry constrained coupling conditions specifying which, if any, free space polarization couples upon perturbation. The degenerated space groups compatible with each unperturbed lattice are exhaustively listed, and in Sec. III.3 the polarization dependence for each mode and each space group is written down by applying Group Theory principles. The result is a catalog of the selection rules for quasi-BICs describing all the unique ways that the six highest symmetry modes of square and hexagonal photonic crystal lattices may be accessed by in-plane symmetry breaking. Finally, in Sec. IV we discuss notable aspects of the catalog and sketch several example device applications. In particular, here and in the accompanying Letter Overvig et al. (2020) we show how the catalog guides device design using successive perturbations to achieve multifunctional control of optical spectra, and introduce a novel class of metasurfaces that use this multifunctional control to spatially control resonant wavefronts.

II Background

II.1 Quasi-Bound States in the Continuum

We begin by reviewing the design process to create a finite-sized DHCG, we explore a BIC in a diatomic lattice artistically depicted in Fig. 1(a). Figure 1(b) and 1(c) define the geometric parameters of the unperturbed and perturbed lattices, respectively. The two “atoms” (here, pillars of Silicon) in the perturbed lattice are identical in height, $H$, and diameter, $D$, and sit in a lattice of period $\sqrt{2}P$, where $P$ is the period of the unperturbed lattice. The perturbation can be quantified as the gaps between atoms: the perturbed gap is $g_{2}=g+\delta$, where $g$ is the unperturbed gap ($g=P-D$) and $\delta$ is the perturbation. The FBZs of the unperturbed and the perturbed lattices are shown in Fig. 1(d) and 1(e), respectively, with high symmetry points defined and the primed coordinates representing the perturbed lattice. The effect of the lattice transformation (taking the period in real space from $P$ to $\sqrt{2}P$ and rotating the basis vectors by $45\degree$) is to shrink the extent of the FBZ and rotate it by $45\degree$. The states belonging to sections of the unperturbed FBZ that lie outside of the new, perturbed FBZ are, by Bloch’s Theorem, equivalent to states within the new FBZ. They are brought into the new FBZ by translation of a reciprocal lattice vector (a process known as Brillouin Zone folding) as depicted graphically in Fig. 1(e) for the shaded area near the $X$ point. The bound modes that were at the $X$ point are now at the $\Gamma$ point (that is, in the continuum) due to the perturbation.

The new modes brought into the continuum may now couple and produce Fano-like sharp spectral features for normally incident light. By construction, the coupling strength is related to the magnitude of the perturbation. It has been shown Overvig et al. (2018) that the coupling strength for small perturbations is of the order of $\delta$. Since the Q-factor of a sharp resonance is inversely proportional to the square of the coupling strength Fan et al. (2003), a symmetry-protected BIC has a Q-factor governed by Overvig et al. (2018); Koshelev et al. (2018)

[TABLE]

where the constant $C$ can vary depending on the mode, geometry, materials, and polarization. Figure 1(f) and 1(g) show the mode profiles for the fundamental mode depicted in Fig. 1(a). Figure 1(h) shows full-wave simulations of the Q-factor of the fundamental mode as a function of perturbation strength, agreeing well with Eq. (1) with $C\approx 6.5P^{2}$. Figure 1(i) contains the band diagram for the perturbed structure calculated by the planewave expansion method (PWEM) using the supercell method, with high symmetry points defined relative to both the unperturbed and perturbed lattices. The modes are calculated in the unperturbed structure following the red arrows in Fig. 1(d), and then artificially folded into the FBZ.

The band structure of the perturbed PCS in Fig. 1(i) can help predict the accuracy of Eq. (1) for finite devices. In an infinite device, a planewave corresponds to a single state (for instance, a mode at the $\Gamma$ point) and the band curvature is irrelevant. However, a finite device excited by a Gaussian beam will behave as some combination of responses excited by the planewaves composing that Gaussian beam. One simple model for predicting the behavior of a finite device is to perform a weighted sum of the spectra corresponding to the constitutive planewaves Overvig et al. (2018). We model a band by a Taylor expansion about the $\Gamma$ point, $\omega_{res}(k)=\omega_{0}+bk^{2}$, where $\omega_{0}$ is the angular frequency of the mode at $k=0$ and $b=\frac{1}{2}\left.\frac{\partial^{2}\omega_{res}}{\partial k^{2}}\right|_{k=0}$. A Gaussian beam with a characteristic spread in wavevector of $\Delta k$ will excite a characteristic spread of frequencies $\Delta\omega=b\Delta k^{2}$. It is natural to expect that if this spread of frequencies is larger than the linewidth of the resonance, $d\omega$, excited in an infinite device by a planewave, the spectral feature will be washed out, lowering the observed $Q$ and invalidating Eq. (1). This suggests a constraint

[TABLE]

In other words, there is an upper limit on the Q-factor attainable in a finite device due to the band curvature near the $\Gamma$ point.

While this simple model does not account for all of the possible finite size effects (e.g., edge effects and a more complex modal structure), the derived constraint suggests that optimizing the band flatness will tend to allow for the most compact devices. In particular, the factor to minimize is $F=\left\lvert b\right\rvert/\omega_{0}$, which serves as a figure of merit when designing a device by computing its band structure. Figure 1(j) maps $F$ as calculated for a variety of diameters and heights (relative to the period) of Silicon pillars sitting on a Silicon Dioxide substrate. While the curvature is different along $\Gamma-M^{\prime}$, or the $k_{x}$ direction, compared to along $\Gamma-X^{\prime}$, or the $k_{a}$ direction, this band is limited by its curvature in the $k_{x}$ direction; we therefore restrict the calculation of $F$ to the band along the $k_{x}$ direction.

We choose a design with the smallest $F$ according to Fig. 1(j) and scale its geometrical parameters by a factor $\lambda/\lambda_{res}$ such that the operating resonant wavelength is $\lambda=1.58\mu m$ for a calculated resonant wavelength, $\lambda_{res}$. Figure 1(k) shows transmittance spectra calculated by full-wave simulations of an infinitely periodic device excited by a planewave of either $x$ or $y$ polarization, demonstrating that this coupling only occurs for $y$ polarization. Figure 1(k) also shows a transmittance spectrum of a device of finite size ($30\mu m\times 30\mu m$) excited by a Gaussian beam with a waist radius of $w_{0}=7\mu m$ calculated by full-wave simulations. The spectral feature remains intact, confirming that the flat band in Fig. 1(i) determined through the optimization shown in Fig. 1(j) allows for compact devices with moderately high $Q\approx 10^{3}$.

Figure 1 overviews the design process of a compact optical device (a two-dimensional DHCG) supporting a sharp spectral feature due to a quasi-BIC. However, this process represented just one high symmetry mode, and explored the behavior as a result of only one specific perturbation. This behavior was shown to be weakly dependent on incident angle (Fig. 1(i)), but strongly dependent on incident polarization (Fig. 1(k)). The key result of this paper is a theoretical description allowing prediction of the polarization behavior (selection rules) of all high symmetry modes and perturbations. With this result, which we call the “catalog of selection rules”, the design process overviewed in Fig. 1 can be summarized as having three steps: (1) Choose a high symmetry mode for its real-space properties (e.g., for its field overlap with the high index material); (2) Optimize the band structure by tuning the parameters in unperturbed structure (as in Fig. 1(i)); (3) Choose a proper perturbation according to the desired selection rules (e.g., targeting $y$ polarization). The catalog serves as a comprehensive guide for step (3), clarifying the wealth of options in conjunction with the choice in step (1) of the desired high symmetry mode; it thereby provides a high-level roadmap for this three step design process resulting in a PCS that confines light in both space and time. This design scheme may be further coupled with computational inverse design techniques Molesky et al. (2018) to reduce the dimensions of the design parameter space to be explored.

II.2 Classification of High Symmetry Modes

The first step to determining the selection rules of perturbed PCSs is to classify the modes present. Since the selection rules arise from symmetry breaking, a mode classification scheme employing the symmetries of the allowed modes is the natural choice. Although the final devices of interest are three-dimensional in nature (having a finite thickness in the out-of-plane direction, $z$), it considerably simplifies the analysis to begin with Maxwell’s equations in two dimensions. In this case, Maxwell’s curl equations decouple into two separate sets of three equations, each set defining modes characterized by either the out-of-plane magnetic field, $H_{z}$, (referred to as TE modes) or the out-of-plane electric field, $E_{z}$ (referred to as TM modes). Each mode is then definable by this single field component. We therefore select, review, and carry out a Group Theory approach detailed in Ref. Sakoda (2005) to classify the modes by in-plane symmetries of the out-of-plane field component. We note that this Group Theory analysis is valid for any materials system, for instance, an array of Silicon pillars, holes in a Silicon slab, or even a metallic structure. For convenience, and comparison to conventional metasurfaces, we first consider arrays of Silicon pillars. But the resulting selection rules are immediately transferable to any other materials system.

Since the fields exist in a periodic lattice, they are characterized by planewaves with magnitudes and directions corresponding to high symmetry points of the reciprocal lattice. When the index contrast is low, this characterization is excellent; when the index contrast is large (such as a Silicon and air system), significant deviations in resonant frequencies occur relative to the low index contrast systems, but the symmetries of the possible modes remain unchanged. The modes can therefore be studied with reference to the extended zone scheme.

Figure 2(a) depicts the extended zone scheme for a square lattice, with notable high symmetry points marked. In particular, the $\Gamma$ points, $M$ points, and $X$ points are labeled by an index pertaining to their distance from the origin, $\Gamma^{(0)}$. These points have point group symmetries $C_{4v}$, $C_{2v}$, and $C_{4v}$, respectively (see Appendix A for the character tables and other relevant Group Theory tables), and the modes decomposable by planewaves corresponding to these points are describable by these point groups. These three sets of points are the highest symmetry points in the reciprocal lattice, and therefore correspond to the modes of interest in the square lattice. The three analogous sets of points in the hexagonal lattice (Fig. 2(b)) are the $\Gamma$ points, $M$ points, and $K$ points.

A Group Theory approach Sakoda (1995, 2005); Yu et al. (2018) predicts the number and nature of the modes from each set of high symmetry points in the extended zone scheme. Figure 2(c) contains a table summarizing the modes possible at each of the high symmetry points in the square lattice. The degeneracy of a set of high symmetry points, $N$, is also the number of modes corresponding to that set. The magnitude of the wavevectors, $k$, of the planewaves of a set will correspond to the expected eigenfrequencies of the modes (however, as noted above, this correspondence is poor in high-index contrast systems). Lastly, the irreducible representations describe the mode symmetries. That is, modes that “transform like” (share all the symmetries of) each irreducible representation listed in an extended zone will be present at that extended zone. Note that the $E$ irreducible representations are doubly degenerate, and so account for two modes.

Figures 3 and 4 depict the TM modes from the first four extended zones of each high symmetry point in the square and hexagonal lattices, respectively. An analogous set exists for TE modes, identical in symmetries (in $H_{z}$ instead of $E_{z}$) but spatially distorted and differing in eigenfrequency. The modes are organized by the extended zone order (columns) and irreducible representation (rows). Reference to the relevant character tables (Appendix A, Fig. A 1(a)) shows that modes labeled by a given irreducible representation transform the same way as the corresponding row in the character table: a $1$ in a column of this row means the mode will be symmetric under the class of operations of that column; a $-1$ means anti-symmetric; a [math] means not symmetric; and a magnitude of $2$ signifies that the mode is degenerate.

Finally, the out-of-plane property of the modes is characterized by the order, $n$, or number of anti-nodes, per atom of the PCS in the $z$ direction. The inclusion of out-of-plane characteristics captures all the relevant features of the modes within the scope of this paper if the PCS has mirror symmetry about an $xy$ plane. However, two-dimensional PCSs with a substrate are known to exhibit chiral behavior: incident circularly polarized light can behave in a manner depending on the handedness Wu et al. (2014); Eismann et al. (2018). The chiral effects of a substrate and vertical symmetry breaking are beyond the scope of this paper, and represent a fruitful avenue for future research. We restrict ourselves to PCSs composed of vertically extruded 2D lattices, and we find that the presence of a low-index substrate (such as glass) generally has little practical effect of this kind (and so can be ignored).

With the in-plane and out-of-plane features of each mode classified, we are motivated to provide a naming scheme. We call a mode:

[TABLE]

where $\psi$ is TM or TE if the mode is characterized by $E_{z}$ or $H_{z}$, respectively, $L$ signifies the reciprocal lattice point (e.g., $\Gamma$), $S$ is the irreducible representation (e.g., $A_{1}$), $m$ is the extended zone order, and $n$ is the out-of-plane order. For instance, the mode in the $B_{2}$ row and $X^{(1)}$ column in Fig. 3, with a single out-of-plane anti-node per unit cell of the PCS would be called $\text{TM}_{X,B_{2}}^{1,1}$, which is the lowest frequency $E_{z}$ mode in this square lattice. $\text{TE}_{X,B_{2}}^{1,1}$ is the mode explored in Fig. 1.

Importantly, we discuss the relationship of the 2D description of the modes and the modes of a finite-thickness PCS, which we assume throughout this Article is simply extruded (its cross-section is invariant) in the $z$ direction. The modes depicted in Fig. 3 and Fig. 4 may be considered as the modes traveling in the $z$ direction in a semi-infinite 2D photonic crystal. In this case, these eigenmodes are described by folding the eigenmodes of an unpatterned isotropic medium (i.e., planewaves), while in a finite PCS, as discussed in Sakoda (2005), the eigenmodes are described by folding the eigenmodes of the unpatterned slab. While the underlying bases are distinct, the resulting modal symmetries are identical. However, polarization mixing occurs in the PCS as a result of the finite thickness, yielding modes that are quasi-TE and quasi-TM, rather than pure. But due to the vertical extrusion, and for normally incident light, the new polarization components cannot introduce or destroy any symmetries, which are described by the point group; that is, the TM (TE) components of the quasi-TE (quasi-TM) modes contain equivalent symmetry properties to the analogous pure-TE (pure-TM) modes, and therefore have no bearing on the selection rules for normally incident light. We refer to quasi-TE (quasi-TM) modes as simply TE (TM) for this reason, and are free to treat the PCS modes as equivalent to their 2D counterparts for the purposes of studying the selection rules at normal incidence.

The modes of a PCS also may have additional, out-of-plane symmetries compared to their 2D counterparts. If there is a $xy$ mirror-plane (at the center of the PCS), then the point groups used here are not sufficient to describe the modes. For instance, the modes of the hexagonal lattice are described by the dihedral group $D_{6h}=C_{6v}\otimes C_{1h}$, where the $C_{1h}$ group accounts for whether the modes are symmetry or anti-symmetric with respect to that $xy$ mirror plane. This symmetry determines, in part, the “handedness” of the Fano resonance Fan et al. (2003) (that is, whether the reflection peak occurs at a redder or bluer wavelength than the reflection dip). However, it has no bearing on the selection rules at normal incidence, and may therefore be ignored for our purposes. Furthermore, in the most practical scenario a low-index substrate is present, breaking this mirror symmetry and leaving only in-plane symmetries. Since inclusion of the mirror plane needlessly doubles the number of modal labels, we exclude it.

Instead, we use the modal index $n$ to refer the out-of-plane characteristics of the PCS modes. When $n$ is odd (even), the decay symmetry is even (odd), meaning that the mode naming scheme in Eq. 3 contains the relevant information about the handedness of the asymmetric lineshape (note that the decay symmetries hold approximately even in the presence of a low-index substrate). Consequently, there is a 1:many correspondence between the modes of a semi-infinite photonic crystal and a PCS: each mode in Figs. 3 and 4 has identical in-plane symmetries to many modes in the PCS that differ only out-of-plane according to the modal index $n$. Such a relationship between accidental BICs supported by an HCG and the vertically propagating waveguide modes of the corresponding 1D waveguide array is described in Chang-Hasnain and Yang (2012), wherein a round-trip phase condition of the vertically propagating modes predicts the dispersion of the BICs; here, the integer multiple of $2\pi$ picked up upon a round trip is the modal index $n$.

Finally, we note that the $K$ point modes in the hexagonal lattice are more complex than the other five high symmetry modes. In particular, the $K$ point in the unperturbed lattice has point symmetry of $C_{3v}$, as evident in Fig. 2(b). However, there are two identical sets of these $K$ points; the set not pictured in Fig. 2(b) can be obtained by reflection about the $k_{y}$ axis. As demonstrated in Fig. 1, the modes only become BICs once folded to the $\Gamma$ point by an appropriate perturbation (depicted in Fig. 5). For $K$ point modes, the analogous perturbation results in a triatomic lattice and therefore triples the number of modes at the new $\Gamma$ point compared to the unperturbed monatomic lattice. At the perturbed $\Gamma$ point, a set of modes originating from each set of $K$ points will mix in pairs. The symmetries of the mixed modes are described by the direct product group $C_{s}\otimes C_{3v}$, corresponding to the relationship of the two sets of $K$ points. This direct product group is isomorphic to (shares the same character table as) the group $C_{6v}$, allowing the modes to be named in $C_{6v}$. Naming the modes according to $C_{6v}$ is inconsistent with the modes in an unperturbed lattice at a single $K$ point, but consistent with the modes upon folding to the $\Gamma$ point and mixing in the relevant perturbed lattice. Since our goal is to study these modes in the perturbed lattice, defining the modes in $C_{6v}$ is the more fruitful choice.

III Selection Rules

By proper periodic perturbation, any of the six classes of high symmetry modes can be accessed from free space if additional symmetry constraints are satisfied. These symmetry constraints can be treated with a Group Theory approach, and result in a catalog detailing how each high symmetry mode classified above couples to free space under a given planar perturbation. In the following, we identify six lattice types chosen to target the six high symmetry modes (Fig. 5), list all the degenerated space groups compatible with those lattice types (Fig. 6), and then derive the selection rules for every case (exemplified in Figs. 7 and 8). The resulting catalogs (Figs. 9 and 10) can be used as a high-level roadmap in the design of planar photonic devices.

III.1 Target Space Groups

The six types of high symmetry modes described above motivate six types of lattices, each one uniquely targeting one of the six high symmetry mode types. For each of these lattices, an exhaustive list of lattices with lower symmetry attainable by planar perturbation is explored. The symmetry degeneration from higher symmetry to lower symmetry will constrain which polarization, if any, may couple to free space for each high symmetry mode.

The six lattice types, depicted in Fig. 5, are named based on the modes they target and whether they begin with square or hexagonal symmetries. For instance, the $Sq_{\Gamma}$ is a monatomic photonic crystal with a square lattice where the perturbation has periodicity equal to that of the unperturbed lattice. This lattice is labeled by $\Gamma$ because it supports none of the other types of modes of interest supported by the square lattice (that is, $M$ and $X$ modes) in the continuum. Figure 5 (top left) depicts an example real space lattice, First Brillouin Zone (FBZ), and band diagram for the $Sq_{\Gamma}$ lattice. The white region in the band diagram is the region of the continuum of interest, wherein only the 0th diffractive order is allowed. We constrain ourselves to the area near the $\Gamma$ point of the white region, where the symmetry-protected BICs can produce sharp spectral features described above.

The $Sq_{M}$ lattice (top middle of Fig. 5), on the other hand, is a photonic crystal with perturbations with periodicity double that of the unperturbed lattice in a single direction. This period doubling (in the $x$ direction in Fig. 5) halves the extent of the FBZ in the $k_{x}$ direction. The shaded portion outside the new FBZ is then translated into the FBZ by a reciprocal lattice vector. As a result, the $M$ point of the unperturbed lattice overlaps with the $\Gamma$ point, bringing the $M$ point modes into the continuum in an analogous way described in the example in Fig. 1. This Brillouin Zone folding also changes the shape of the 0th order diffraction region of the band diagram. The $\Gamma$ point will now have both the modes at the unperturbed $\Gamma$ point as well as at the unperturbed $M$ points. The $Sq_{M}$ lattice is the only lattice in Fig. 5 to bring the $M$ point modes of a square lattice into the continuum, motivating its name. The remaining lattices target $X$ modes of the square lattice ($Sq_{X}$, which is the lattice type explored in Fig. 1), and $\Gamma$, $M$, and $K$ point modes of the hexagonal lattice ($Hex_{\Gamma}$, $Hex_{M}$, and $Hex_{K}$, respectively) in an analogous way. Notably, two distinct regions are folded into the FBZ of the $Hex_{K}$ lattice: as discussed above, two sets of modes are folded to the $\Gamma$ point, one from each distinct $K$ point.

We note that the six lattices chosen in Fig 5 are not an exhaustive set: lattices with any number of atoms per unit cell are possible. Ordering this list of lattices by number of atoms per unit cell, the six chosen lattices are the lowest order lattices uniquely targeting the six high symmetry modes of interest. Appendix B describes three examples of higher order lattices. The approach described in what follows may be used to generate the catalog for any higher order lattice.

Next, the space groups of degenerated lattices that are compatible (attainable through perturbation) with each lattice are determined and reported in Fig. 6, the space group compatibility table. First, all $17$ “wallpaper groups” are listed and categorized by the compatible lattice family (e.g., “Rectangular”). The point group of each of these is given for reference. Then, for each target lattice type (e.g., $Sq_{M}$) the space groups compatible with the lattice class (i.e., square or hexagonal) are listed. For instance, the space group $p6mm$ requires hexagonal tiling, and is therefore omitted as a possibility for any $Sq$ lattice. Likewise, a $Hex$ lattice cannot be perturbed into a square lattice without distortion of the lattice vectors, but it can be perturbed into a rectangular lattice. The space groups of the square lattice family are therefore omitted from all $Hex$ lattices, but those of the rectangular lattice family are not.

Next, the glide reflection operation (a reflection and a translation by a fraction of a unit cell) is tested for each type of lattice. Glide reflections are present in only some two-dimensional space groups and are not compatible with all of the six target lattice types. As an example, it is quickly found by inspection that the $Sq_{\Gamma}$ lattice does not support glide symmetries in directions other than along the diagonals (more rigorously, in the language of crystallography, monatomic PCSs are incompatible with non-symmorphic space groups). This excludes the space groups $pg$, $pgg$, and $pmg$, which are correspondingly greyed out in the column for $Sq_{\Gamma}$. For the same reason, for the $Sq_{M}$ lattice, glide planes along the direction where the lattice is unperturbed (and therefore monatomic) are incompatible (the $y$ direction in Fig. 5). Additionally, all diagonal glides are incompatible for the $Sq_{M}$ lattice because they correspond to reflection axes that are not included in the point group of any $Sq_{M}$ lattice. This excludes $cm$, $pg$, and $cmm$ for $Sq_{M}$, which are greyed out accordingly. There are no such constraints for the $Sq_{X}$ lattice, which can be degenerated into a lattice of any space group (except the hexagonal ones). Similar arguments can be made for the $Hex$ lattices, and the results are reported in Fig. 6.

Finally, it must be noted that there exist multiple high symmetry points in each real space lattice. These are given names in Fig. 5 for each case. For instance, the $Sq_{\Gamma}$ lattice has two points having the full symmetry of the $C_{4v}$ point group, named $\gamma$ and $\xi$. Both are perfectly acceptable to choose as the reference point: in the mode naming scheme in Sec. II.2, the $\xi$ point is the reference point, but the modes may all be renamed according to the $\gamma$ point if desired. Similarly, a degenerated space group may choose either of these points to have in common with the unperturbed lattice. Generally speaking, every degenerated space group may be tried with each of the high symmetry points in common with the unperturbed lattice, thereby allowing for more than a single unique example of each space group in each lattice. For instance, there are three distinct $cmm$ space groups in the $Sq_{X}$ catalog (see Fig. 9): one with the $\gamma$ point as the high symmetry point in common, one with $\xi$, and the last with $\mu$ (which is the space group of the perturbation in Fig. 1). As shown in Fig. 9, though these have identical space groups, they do not have identical selection rules because they are attained through distinct perturbations. Therefore, to determine all of the unique symmetry degenerations possible, an attempt is made to construct each compatible space group (Fig. 6) with each high symmetry point in common between the unperturbed and perturbed lattice. The successful attempts comprise the set of all degenerated lattices compatible with those chosen in Fig. 5. This proof by exhaustion is omitted here. For each of these degenerated lattices, the modes derived in the previous section can be studied, and their selection rules derived. The results are tabulated in Figs. 9 and 10 following the methods derived in Secs. III.2 and III.3.

III.2 Deriving the Coupling Condition

To derive the selection rules reported in Figs. 9 and 10, we study the end-fire coupling of free-space light normally incident to a semi-infinite 2D photonic crystal supporting the modes depicted in Fig. 3 and Fig. 4. As described in Sec. II.2, there is a 1:many correspondence preserving in-plane symmetries between the modes excited by this end-fire coupling and the modes excited in a finite-height PCS (under the assumptions that the cross-section is invariant along the $z$ direction, whose substrate and superstrate are isotropic media). At normal incidence, the two cases therefore have identical selection rules with respect to in-plane symmetry breaking, and so we may proceed with the simpler case of end-fire coupling. The multiplicity of modes of the PCS differentiated by the modal index $n$ have identical selection rules; it is only the in-plane symmetries that are relevant. In particular, we determine under which conditions the end-fire coupling coefficient, $\gamma_{e}$, is non-vanishing:

[TABLE]

where $\hat{\boldsymbol{z}}$ is the unit vector in the $z$ direction and bold symbols refer to vector quantities. We exclude a normalization factor for simplicity (it does not affect whether $\gamma_{e}$ vanishes) and evaluate the integral over the area, $A$, of a unit cell. Here, we take the incident field to be a normally incident planewave, with electric field

[TABLE]

and magnetic field

[TABLE]

where $\eta_{0}$ is the impedance of free space. And the mode has electric field

[TABLE]

and magnetic field

[TABLE]

where $\beta$ is the propagation constant satisfying the dispersion relation $\omega=c\beta$. Evaluating the cross-products in Eq. (4), the free-space coupling coefficient is written as

[TABLE]

Using Maxwell’s curl equations, the in-plane components ($E_{x},E_{y},H_{x},H_{y}$) are replaced with the out-of-plane components ($E_{z},H_{z}$) to both simplify the equation and allow the previous mode classification scheme (based on the out-of-plane field components) to straightforwardly apply. The resulting form is

[TABLE]

where

[TABLE]

and

[TABLE]

with the replacement $\varepsilon_{r}(x,y)=\varepsilon(x,y)/\varepsilon_{0}$ as the relative permittivity. This can be written more compactly as

[TABLE]

where angled brackets indicate integration over a unit cell. While it is possible to proceed with this form by considering the symmetries of each component, it is considerably simpler and more informative to reduce this to individual choices of incident polarization (e.g., choose $A_{x}=0$) and mode type (i.e., choose either TM modes or TE modes). In this case, we write

[TABLE]

where $\psi$ is TM or TE, $\partial_{i}$ refers to the partial derivative in a relevant high-symmetry direction ($i=x,y,a,b$), and $c_{j}$ is $c_{1}$ when $\psi$ is TM and $c_{2}$ when $\psi$ is TE.

III.3 Determining the Selection Rules

We now apply Eq. (14) to the modes supported by a 2D photonic crystal. The modes in the unperturbed lattice can be described as the eigenvectors of an eigenvalue equation

[TABLE]

where the superscript marks reference to the unperturbed eigenvalue problem. We are interested in particular in the $\psi^{0}$ that are uncoupled to free space (i.e., $\psi^{0}$ for which the integral in Eq. (14) vanishes) due to symmetry. To proceed we apply perturbation theory to determine any non-vanishing terms present in the generalized eigenvalue problem of a degenerated lattice:

[TABLE]

where $H=H^{0}+V$ is perturbed by the perturbation operator, $V$, and $\psi=\psi^{0}+\psi^{1}$ is the perturbed field profile with $\psi^{1}$ as the first order correction. First order perturbation theory gives the form of the $n$th mode $\psi_{n}^{1}$ as

[TABLE]

That is, the perturbed portion of the field is a superposition of the unperturbed fields. (Note that Eq. (17) is the non-degenerate form of perturbation theory, but it can be applied to degenerate states as well if the “correct” orthogonal linear combination of states is known ahead of time. Since these will correspond to a high symmetry direction of the perturbed lattice, they are easy to predict; we therefore use this form to apply to degenerate modes.) However, inspection of the coefficient, $\left\langle(\psi_{m}^{0})^{*}V\psi_{n}^{0}\right\rangle$, in front of each $\psi_{m}^{0}$ reveals that not all $\psi_{m}^{0}$ will contribute: many will vanish due to symmetry.

The process of determining $\psi_{m}^{0}$ that contribute can be clarified and expedited in the language of Group Theory. Specifically, if we can determine the irreducible representations of each factor within the integrand, we can find the symmetries of the total integrand by computing the direct product of those irreducible representations. A direct product is an abstract way to obtain the symmetries of the product of two functions $f$ and $g$: if $h=fg$, the symmetries of $h$ may be obtained by performing the direct product of the representations of $f$ and $g$. That is, $\Gamma_{h}=\Gamma_{f}\otimes\Gamma_{g}$, where $\Gamma_{h}$ is the representation of $h$ in some point group. Since the fields have been classified already in terms of their irreducible representations, we write the irreducible representation of $\psi_{m,n}^{0}$ as $\Gamma_{m,n}$ and simply refer to Figs. 3 or 4. Then, we must determine the irreducible representation of $V$, or $\Gamma_{V}$, which can be achieved following a process described below. We finally write the direct product as $\Gamma_{integrand}=\Gamma_{m}\otimes\Gamma_{V}\otimes\Gamma_{n}$.

A necessary condition for this integral to be non-vanishing is that this direct product must contain a component that transforms as a constant: the sinusoidal components do not contribute upon integration over a unit cell. Since a constant is fully symmetric (that is, it transforms as $\Gamma_{1}$, which is the highest symmetry irreducible representation in every point group), this condition is identical to saying that $\Gamma_{integrand}$ must contain $\Gamma_{1}$. Note that this condition is necessary, but not sufficient. For instance, a cosine transforms as $\Gamma_{1}$ about the origin, but integrates to zero over a period. We can therefore say that $\left\langle(\psi_{m}^{0})^{*}V\psi_{n}^{0}\right\rangle$ is non-vanishing only if $\Gamma_{m}\otimes\Gamma_{V}\otimes\Gamma_{n}=\Gamma_{1}+...$.

A direct product is easily calculated by referring to the direct product table of the relevant point group (see Appendix A, Fig. A 2). A notable feature of these tables is that two irreducible representations $\Gamma_{i}$ and $\Gamma_{j}$ satisfying $\Gamma_{i}\otimes\Gamma_{j}=\Gamma_{1}+...$, also satisfy $\Gamma_{i}=\Gamma_{j}$. Consequently, we can reframe the condition on the integrand, $\Gamma_{m}\otimes\Gamma_{V}\otimes\Gamma_{n}=\Gamma_{1}+...$ to be $\Gamma_{m}=\Gamma_{V}\otimes\Gamma_{n}$. In other words, a field $\psi_{m}^{0}$ contributes to the perturbed field $\psi_{n}^{1}$ only if $\Gamma_{m}=\Gamma_{V}\otimes\Gamma_{n}$. Since $\psi_{n}^{1}$ will transform as the components comprising it (that is, the $\psi_{m}^{0}$ with non-vanishing integrals), we finally arrive at the conclusion that $\Gamma_{\psi_{n}^{1}}=\Gamma_{V}\otimes\Gamma_{n}$. Since the index $n$ refers to any particular mode of interest, we may drop it:

[TABLE]

That is, the first order perturbed field profile transforms as the direct product of the irreducible representations of the perturbation operator and the unperturbed field profile in question.

Understanding the symmetries of the perturbed portion of the wavefunction allows us to simplify the free-space coupling condition, Eq. (14):

[TABLE]

Since $c_{j}$ is a function of in-plane permittivity distribution of a perturbed lattice, it is natural to expect that it has a portion that transforms like $H^{0}$, which we call $c_{H^{0}}$, and a portion that transforms like $V$, which we call $c_{V}$. We can write these portions explicitly to first order using the binomial approximation. Taking $c_{2}=c_{H^{0}}+c_{V}$, for instance, the unperturbed portion is written

[TABLE]

and the perturbed portion is written

[TABLE]

where $\varepsilon_{r}(x,y)=\varepsilon_{r}^{0}(x,y)+\delta\varepsilon_{r}(x,y)$ is decomposed into the unperturbed portion, $\varepsilon_{r}^{0}(x,y)$, and perturbed portion, $\delta\varepsilon_{r}(x,y)$, of the relative permittivity. Since $\varepsilon_{r}^{0}(x,y)$ transforms as $\Gamma_{1}$ by construction of the unperturbed lattices, it is evident that $c_{H^{0}}$ transforms as $\Gamma_{1}$: for an even function $f(x)$, the function $1/(1-f(x))$ is also even. Decomposing the factors in $c_{V}$, it is clear that it transforms as $\Gamma_{V}$ because $\delta\varepsilon_{r}(x,y)$ transforms as $\Gamma_{V}$ by definition of the perturbation, and the remaining factors in $c_{V}$ transform as $\Gamma_{1}$ (which acts as the identity in direct products). A similar argument reveals the equivalent result for $c_{1}$.

We therefore write $\gamma_{e}$ as the sum of four terms:

[TABLE]

The first term vanishes for symmetry-protected BICs. As described above, $\Gamma_{H^{0}}=\Gamma_{1}$, and so, using Eq. (18), the second term’s integrand transforms as $\Gamma_{1}\otimes\Gamma_{\partial_{i}}\otimes\Gamma_{V}\otimes\Gamma_{\psi^{0}}=\Gamma_{\partial_{i}}\otimes\Gamma_{V}\otimes\Gamma_{\psi^{0}}$. The third term’s integrand straightforwardly transforms as $\Gamma_{V}\otimes\Gamma_{\partial_{i}}\otimes\Gamma_{\psi^{0}}$, identical to the second term (inspection of Appendix A, Fig. A 2 shows that the direct products in question commute). The fourth term vanishes to first-order, because it is the product of two factors of the perturbation. We are therefore left with two generally non-vanishing terms whose integrands transform identically.

As before, $\gamma_{e}$ is non-vanishing only if the integrand has a component that transforms like $\Gamma_{1}$. We therefore arrive at the symmetry constrained coupling condition:

[TABLE]

Since a partial derivative in the $i$ direction transforms like a vector in that direction, it also transforms the same as a free space polarization $i$. The physical interpretation of the coupling condition, Eq. (23), then, is that the symmetries of the perturbed part of the field (i.e., $\Gamma_{V}\otimes\Gamma_{\psi^{0}}$) must match the symmetries of a free space polarization (i.e., $\Gamma_{\partial_{i}}$). That is to say, the perturbed field must have a net dipole moment to couple to free space.

The coupling condition is equivalent to considering whether the integral

[TABLE]

vanishes. This form justifies a convenient and insightful graphical method Kilic et al. (2008) of determining whether $\gamma_{e}$ is non-zero, without directly determining $\psi^{1}$, which is not obvious at first glance at Eq. (17). The perturbed mode can be simply drawn by altering the magnitude of the unperturbed field according to the shape and sign of the perturbation. Then, this new perturbed field is decomposed into the unperturbed portion and the perturbed portion (corresponding to $V\psi^{0}$). Taking the derivative amounts to treating the product $V\psi^{0}$ as “charges” and the gradient as the “moment”; then, if there is a net dipole moment, the mode couples to the corresponding free space polarization. Figure 7 depicts this process for determining the selection rules of the $\text{TM}_{X,S}^{m,n}$ modes in a $Sq_{X}$ lattice with a $cmm$ space group (the same used in Fig. 1). The polarization depicted corresponds to the out-of-plane field component. That is, if $\psi^{0}$ is a TE (TM) mode, the polarization depicted describes the magnetic (electric) polarization of free space that couples to $\psi$. Figure 7 therefore correctly predicts the polarization dependence seen in Fig. 1 for $\text{TE}_{X,B_{2}}^{1,1}$.

A more expedient method to generate the selection rules, however, is to determine the irreducible representations present in $V$ and then employ the direct product tables (see Appendix A, Fig. A 2(a) and A 2(b)) to immediately write the selection rules for all modes present at the $\Gamma$ point of the perturbed lattice. This is done by (1) finding the point group in common among $V$ and $\psi$, (2) writing the irreducible representations of each factor in that point group, and then (3) determining if the direct product $\Gamma_{V}\otimes\Gamma_{\psi}$ matches the irreducible representation of a free space polarization (which are reported for each relevant point group in Appendix A, Fig. A 2(c)).

The irreducible representations of $V$ can be found by conventional Group Theory methods if required, but are generally apparent by inspection. Figure 8 depicts the decomposition of $V$ for the same space group as Fig. 7. The process is simplified by properly choosing $H^{0}$ such that $V$ transforms as simply as possible. For instance, $H^{0}$ is written as a circle with permittivity $\varepsilon$ shadowing a square cross oriented in the $x,y$ directions with permittivity $\varepsilon/2$, as shown in Fig. 8(a). It is then clear to see that the $V$ depicted obtains $H$ upon addition of $H^{0}$.

Next, $V$ can be decomposed into two portions, one (called $V_{\Gamma}$) with the periodicity of the unperturbed lattice, and one (called $V_{X}$, here) with the periodicity of the perturbed lattice (Fig. 8(b)). Importantly, Eq. 23 refers only to functions as they exist in the unperturbed FBZ, in which modes characterized by the $X$ point are orthogonal to modes characterized by the $\Gamma$ point. Consequently, $V_{\Gamma}$ contributes only to the modes at the $\Gamma$ point in the unperturbed lattice, while $V_{X}$ contributes only to the modes at the $X$ point in the unperturbed lattice. The point group of both $V_{\Gamma}$ and $V_{X}$ is the same as the point group of the space group, $C_{2v}^{v}$. Referring to the character table of $C_{2v}^{v}$ (Fig. 8(c) or Appendix A, Fig. A 1(a)), it is readily apparent that $V_{\Gamma}$ transforms as $A_{1}$ and $V_{X}$ transforms as $B_{2}$.

Finally, the coupling constraint (Eq. (23)) is evaluated. However, since the modes of interest are defined in a higher group than that of $V$, we must first determine how they degenerate into the lower group. This can be done by referring to the symmetry degeneration tables (see Appendix A, Fig. A 1(c)). Then, the direct products $\Gamma_{V_{\Gamma}}\otimes\Gamma_{\psi}$ are taken with reference to the direct product table for the point group $C_{2v}^{v}$. Since $x,y$ polarized planewaves transforms as $B_{1},B_{2}$ in $C_{2v}^{v}$, the modes for which the product $\Gamma_{V}\otimes\Gamma_{\psi}=B_{1},B_{2}$ couple to $x,y$ polarization, respectively. Figures 8(d) and 8(e) show a worksheet of this process. It bears repeating that this polarization corresponds to that of the out-of-plane field component. For example, if $\psi$ is a $H_{z}$ mode, then $x$ refers to the $H_{x}$ component of the free space planewave, corresponding to $y$ polarized light as conventionally defined by the electric field. The resulting selection rules are in agreement those derived using the graphical method in Fig. 7, but derivation required a single diagram to decompose $V$ instead of one for each mode, and straightforwardly gave the selection rules for the $\Gamma$ point modes as well (for which the graphical method would require another set of diagrams).

The method detailed throughout this section may be summarized as follows. First, the unique space groups compatible with each lattice type are determined by exhaustion. Second, for each of these space groups, the perturbation is split into $V_{\Gamma}$ and $V_{L}$, where $L$ is the high symmetry point of the unperturbed reciprocal lattice that is folded to the $\Gamma$ point. Third, the irreducible representation of each portion of the perturbation is determined. Fourth, the coupling condition, Eq. (23), is evaluated for each mode, using $V_{\Gamma}$ for $\psi_{\Gamma,S}^{m,n}$ and $V_{L}$ for $\psi_{L,S}^{m,n}$; the matching polarization (if any) is marked down by reference to Fig. A 2(c).

With this method, a catalog for each of the six lattice types is generated, for each of the compatible degenerated space groups described above. The catalogs are given in Fig. 9 for square lattices, and in Fig. 10 for hexagonal lattices. An entry of the catalog lists the space group by name, the point group used to describe $V$, the irreducible representations of $V_{\Gamma}$ and $V_{L}$, an example visualization of the degenerated lattice (using the “keyhole” motif Kilic et al. (2008)), and the selection rules for all the modes present at the $\Gamma$ point in the perturbed lattice.

Note that the selection rules for the two-fold cyclic space group, $p2$, in the catalog are specified by some angle, $\phi$ or $\theta$, which are ill-defined relative to the lattices’ axes: the polarization angle must be numerically determined, will generally change with the magnitude of the perturbation, and may differ between TM and TE modes of the same symmetry. However, for small perturbations, two modes controlled by the same $p2$ perturbation and specified by $\phi$ will be excited by the same polarization angle $\phi$; $\theta$ denotes the angle orthogonal to $\phi$. Note that the Group Theory approach in Fig. 8 can only say that some polarization couples to the mode, but cannot specify $\phi$ and $\theta$; for this, the diagrammatic approach in Fig. 7 is used. The selection rules for $p1$ are ill-defined in a similar way: some polarization couples with a direction unconstrained by group theory, and are therefore specified as $any$. The remaining cyclic space groups, $p3$, $p4$, and $p6$, only allow access to degenerate modes in a polarization independent manner, and so for simplicity are specified as $x,y$.

IV Discussions and Applications

The process described above lays out the derivation of the selection rules for two-dimensional PCSs with in-plane perturbations applied. The resulting catalog, split into Figs. 9 and 10, contains a great amount of information and warrants further discussion and exploration. In particular, a few unique features present in the catalog readily motivate device applications not possible in the simpler one-dimensional PCSs.

For instance, due to the two-dimensional nature of the device, the band structure can be optimized in both in-plane directions, allowing for full optimization of the band structure and thereby optimally compact devices. One-dimensional structures (e.g., devices composed of rectangular grating fingers, invariant in one in-plane direction) can be understood as a special case of a subsection of the $Sq_{M}$ lattice, but with limited to no control over the behavior along the direction of the grating fingers. Additionally, the higher in-plane symmetry of two-dimensional structures means the presence of degenerate $E$-type modes (“partner” modes with identical eigenfrequencies that couple to orthogonal polarizations), which do not exist in one-dimensional structures. This allows for compact, polarization independent devices such as filters and modulators to be designed. The manipulation of degenerate modes is therefore of considerable technological interest. Last, we note a parent-child relationship between higher and lower order space groups within the catalog, and find that child space groups constructed by successively adding distinct parent space groups result in optical control with independent degrees of freedom introduced by the parent space groups. We show that this principle enables controlling a large number of parameters characterizing an optical spectrum, well surpassing the state-of-the-art.

With these considerations in mind, Sec. IV.1 discusses the degenerate modes that exist in monatomic and multi-atomic PCSs, providing a comprehensive set of options for polarization-independent devices using the catalog. Section IV.2 details a device application motivated by this discussion, demonstrating that the degenerate fundamental modes $Hex_{K}$, as controlled by three successive perturbations, are suitable for Terahertz generation via four-wave mixing. As a second application of successive perturbations apparent from studying the entries of the catalog, Sec. IV.3 shows the potential for a PCS on a stretchable substrate to enable mechanically tunable optical lifetimes. Finally, Sec. IV.4 reports the discovery of a geometric phase associated with circularly polarized light coupling into a $p2$ space group of the $Sq_{M}$ lattice, controllable simultaneously with the Q-factor of the resonance. We then show how spatially varying this geometric phase enables a novel class of photonic devices in which the outgoing Fano resonant wavefront is spatially tailored while the non-resonant light is left unaffected.

IV.1 Degenerate Modes

We first consider the nature of degenerate modes. The degenerate modes generally transform as partners of a degenerate irreducible representation (e.g., $E$ in $C_{4v}$), which are written $E^{x}$ and $E^{y}$ corresponding to their dipole moments. Because of this dipole moment, the $E^{x}$ and $E^{y}$ modes in the unperturbed $Sq_{\Gamma}$ lattice generally couple to free space (and the $E_{1}$ modes couple to free space in the unperturbed $Hex_{\Gamma}$ lattice). In other words, the integral

[TABLE]

which has an integrand that transforms as $\Gamma_{1}\otimes\Gamma_{\partial_{i}}\otimes\Gamma_{\psi^{0}}$, is non-vanishing for $E$ modes because $\Gamma_{\partial_{i}}=E$ in $C_{4v}$ (likewise, $\Gamma_{\partial_{i}}=E_{1}$ in $C_{6v}$). The coupling can numerically vanish for certain combinations of angle, polarization, and optical materials, but since these are reasons unrelated to the symmetry arguments above, they are accidental BICs; all modes other than the $\psi_{\Gamma,E}^{m,n}$ for $Sq_{\Gamma}$ and $\psi_{\Gamma,E_{1}}^{m,n}$ for $Hex_{\Gamma}$ in the mode classification are symmetry-protected BICs. For the modes already coupled to free space in the unperturbed lattice, the only significant impact a perturbation has is to split the degeneracy upon symmetry degeneration (for instance, perturbing a lattice with $C_{4v}$ down to $C_{2v}$). In other words, if the lattice is made structurally birefringent, the $E^{x}$ and $E^{y}$ modes will degenerate into irreducible representations in a lower order point group with different eigenfrequencies, but the coupling rate to free space will generally be changed to a negligible degree. For this reason, Figs. 9 and 10 simply label the corresponding entries $x,y$.

Of more interest here are the degenerate symmetry-protected BICs. The $M$ point modes of either square or hexagonal lattices have no such degenerate modes because the $C_{2v}$ point group has no degenerate irreducible representation (a rectangle is not identical in the $x,y$ directions). However, the $\psi_{\Gamma,E_{2}}^{m,n}$ for $Hex_{\Gamma}$, $\psi_{X,E}^{m,n}$ for $Sq_{X}$, and $\psi_{K,E_{1}}^{m,n}$ and $\psi_{K,E_{2}}^{m,n}$ for $Hex_{K}$ are degenerate symmetry-protected BICs. Therefore, a polarization insensitive filter or modulator must use one of these lattices in order to utilize the advantages of quasi-BICs (that is, a Q-factor controllable by Eq. (1) independent of the band structure). We consider the degenerate modes in each of these three lattices in turn.

The $Hex_{\Gamma}$ supports the $E_{2}$ modes, which are uncoupled to free space in the absence of a perturbation ($E_{1}\otimes E_{2}=B_{1}+B_{2}+E_{1}$, which does not contain $\Gamma_{1}$ in $C_{6v}$) and are degenerate BICs. Reference to the catalog shows that reducing the symmetry to $C_{3v}$ or lower may allow free-space coupling to these modes. A polarization independent filter or modulator with Q-factor following Eq. (1) could be made utilizing the $E_{2}$ modes of a $Hex_{\Gamma}$ lattice according to either the $p31m$ or $p3m1$ entry of the catalog.

The $Sq_{X}$ lattice supports degenerate modes that are bound in the unperturbed lattice. Upon perturbation, they are brought to the $\Gamma$ point, allowing coupling to free space at normal incidence. Several space groups in the $Sq_{X}$ lattice leave these modes uncoupled in the continuum, making them BICs, while most others allow coupling, making them quasi-BICs. The space groups with $C_{4v}$ and $C_{4}$ leave the eigenfrequencies degenerate, while lower order symmetry groups introduce birefringent behavior. Therefore, a polarization independent filter or modulator with Q-factor following Eq. (1) may be made utilizing the $E$ modes of a $Sq_{X}$ lattice according to any of the $p4m$, $p4g$, or $p4$ entries of the catalog as reference.

An interesting feature of the catalog is the prediction of coupling of the $E^{x}$ partner of the $E$ modes of a $Sq_{X}$ lattice to either $x$ or $y$ polarized light (equivalently, $y$ polarized light may couple into either the $E^{x}$ or $E^{y}$ partner). Compare, for instance, the $p4m$ with $V_{X}=A_{1}$ (Fig. 11(a)) and $p4g$ with $V_{X}=A_{2}$ (Fig. 11(b)) space groups in the $Sq_{X}$ lattice. Figures 11(c) and 11(d) depict, for the $p4m$ and $p4g$ cases respectively, the field profiles calculated by full-wave simulations at the frequency of the $\text{TE}_{X,E}^{1,1}$ modes excited by $y$ polarized light (magnetically $x$ polarized light). The former shows that the magnetically $x$ polarized planewave couples to the $E$ mode with the apparent dipole in the $y$ direction (that is, the $E^{y}$ partner, as defined in Fig. 7), while for the latter it couples to the $E$ mode with the apparent dipole in the $x$ direction (that is, the $E^{x}$ partner). This mode “twisting” is written in the catalog by writing $y,x$ for $p4m$, in contrast to the entry of $x,y$ for $p4g$, and it is easily predicted by tracking how the partners degenerate to $C_{2v}$ in the worksheet of Fig. 8, or by using the diagrammatic approach in Fig. 7. This is a phenomenon that does not occur in the more often studied $E$ or $E_{1}$ degenerate modes of the $Sq_{\Gamma}$ or $Hex_{\Gamma}$ lattices, and is thus representative of the larger range of behaviors present in multi-atomic lattices. Especially notable is that the dependence of the polarization angle of the incident light on the in-plane orientation of elliptical structures suggests that a geometric phase is associated with this coupling. This will be explored in Section IV.4.

Next, the $Hex_{K}$ lattice is a special case, having degenerate modes of two distinct types. The unit cell is a trimer, and therefore contains three times the number of modes as the unperturbed lattice. As described in Sec. II.2, two sets of modes (corresponding to $K_{1}$ and $K_{2}$ points in the FBZ in Fig. 5) are folded to the $\Gamma$ point. Because of the symmetry of the unperturbed lattice, a mode in one of the sets has a counterpart in the other set with equal eigenfrequency. These pairs mix at the perturbed $\Gamma$ point, producing a final set of modes describable in $C_{6v}$. The lattice therefore supports $E_{1}$ and $E_{2}$ modes newly brought to the $\Gamma$ point by the perturbation, analogous to the $E_{1}$ and $E_{2}$ modes in the $Hex_{\Gamma}$ lattice.

However, inspection of the band diagram for the $Hex_{K}$ lattice near the $\Gamma$ point in Fig. 5 (bottom right) shows many more degeneracies than explainable by the $E_{1}$ and $E_{2}$ modes alone. All of the newly folded modes, in fact, are degenerate, despite the mode classification scheme in Fig. 2 predicting the presence of modes not describable by $E$ irreducible representations. In particular, the fundamental modes are degenerate (last band diagram in Fig. 5), but have irreducible representations $A_{1}$ and $B_{1}$. Although visibly quite different (see the $K^{(1)}$ column of Fig. 4), and having distinct symmetries in $C_{6v}$, they are nonetheless identical in eigenfrequency. This degeneracy is born of the trimerization of the lattice: a pair of modes with the same eigenfrequency are superposed upon translation to the $\Gamma$ point, and can be superposed either in phase or out of phase, producing a pair of distinct modes with identical eigenfrequency.

One consequence of this is that in the $p31m$ $Hex_{K}$ lattice with point group $C_{3v}^{d}$, the $A_{1}$ and $B_{1}$ modes form a degenerate pair that together correspond to a spectral feature that is polarization insensitive (a similar behavior is seen in the $p3$ lattice where $V_{K}$ transforms as a partner of the $E$ irreducible representation: the fourth entry in the $Hex_{K}$ lattice). This reveals another way to consider this degeneracy: the $A_{1}$ and $B_{1}$ modes are partners of the $E$ irreducible representation of the $C_{3v}^{d}$ point group, defined about the $\kappa$ point in the real space lattice (as defined in the $Hex_{\Gamma}$ lattice in Fig.  5). Because the $\gamma$ point has the full symmetries of $C_{6v}$, this description of the modes about the $\kappa$ point misses relevant symmetries; nevertheless, the ability to describe them in $C_{3v}^{d}$ as partners of the $E$ modes means their eigenfrequencies must be identical. Upon further symmetry degeneration (for instance, to $cmm$, or $cm$ in the $Hex_{K}$ lattice), the modes behave differently, splitting in both eigenfrequency and polarization dependence. This behavior is unique to the $Hex_{K}$ lattice in Figs. 9 and 10 because it is the only lattice with more than two atoms per unit cell. Higher order lattices, such as those shown in Fig. B 1, may exhibit similar behavior.

IV.2 Application: Terahertz Generation

We explore the degenerate fundamental modes of the $Hex_{K}$ lattice to aid in generating Terahertz (THz) frequencies through nonlinear processes enhanced by optical resonances. Sketching the design of such a device is a useful exercise to demonstrate the utility and an example use of the catalog. Figure 12(a) depicts a schematic of the device, with a $Hex_{K}$ lattice made of Silicon pillars in the gap of a bowtie antenna resonant to a THz frequency. Figure 12(b) shows an example spectrum of the PCS portion of this device, showing two closely spaced resonances at $\lambda_{1}=3.147\mu m$ and $\lambda_{2}=3.161\mu m$, both excited by $y$ polarized light and associated with the split degeneracy of the $\text{TE}_{K,A_{1}}^{1,1}$ and $\text{TE}_{K,B_{1}}^{1,1}$ modes. If optical power is normally incident at pump wavelengths $\lambda_{a}$ and $\lambda_{b}$ such that $\lambda_{a}=\lambda_{1}$ and $\lambda_{b}=\lambda_{2}$, and a low-frequency bias (corresponding to a radiofrequency with wavelength $\lambda_{3}$) is electrically applied across the antenna, four-wave mixing will produce photons at a THz wavelength with improved efficiency (compared to a bulk material) due to the enhanced light-matter interactions from the PCS and antenna resonances. The case shown in Fig. 12(b) corresponds to $\lambda_{4}=711\mu m$, but Fig.  12(c) shows that $\lambda_{4}$ can be easily tuned by the radius of the central pillar, $R_{1}$. Figure 12(d) confirms that the figure of merit (defined below) is indeed maximum at $\lambda_{4}$ when the pump photons have wavelengths of $\lambda_{1},\lambda_{2}$ ($1/\lambda_{3}=0$ for simplicity, here; it may generally be used to finely and actively tune $\lambda_{4}$).

A key advantage of using these degenerate modes is the unique robustness of the control of both the spectral spacing and linewidths of the resonances. Since the modes are degenerate in the unperturbed lattice, they are necessarily closely spaced in a weakly perturbed lattice. Then, by controlling the radius of one of the pillars, the frequency spacing can be finely tuned. The spectral map in Fig.  12(c) shows the impact of tuning the radius of the central pillar, $R_{1}$, depicting a classic anti-crossing behavior Fan et al. (2003) as the resonance spacing changes. This utilization of degenerate modes offers a considerably more robust control of closely spaced resonances compared to relying on controlling two unrelated resonances by tuning geometric parameters: fine tuning the separation of two unrelated resonances is highly sensitive to fabrication errors, while the split degeneracy here is guaranteed by symmetry.

Notably, $\lambda_{4}$ (or the spacing of the resonances) can be tuned largely independently of the linewidths of the resonances. This is easily understood by considering the portion of the perturbation, $V_{1}$, that corresponds to changing $R_{1}$. Depicted in Fig. 12(e), a lattice perturbed by $V_{1}$ alone produces a $Hex_{K}$ lattice with the $p6mm$ space group; reference to the catalog (Fig.  10) reveals that no coupling to the target modes is introduced by this perturbation. Tuning the radius of the central pillar therefore does not affect the coupling to first order. Then, the addition of $V_{2}$ degenerates the space group to $cmm$, which couples the $\text{TE}_{K,B_{1}}^{1,1}$ mode to the magnetic $x$ polarization, but not the $\text{TE}_{K,A_{1}}^{1,1}$ mode. Finally, the addition of $V_{3}$ creates a lattice with the $cm$ space group with the $C_{s}^{v}$ point group, allowing coupling the $\text{TE}_{K,A_{1}}^{1,1}$ mode to the magnetic $x$ polarization. Notably, if the other $cm$ space group (with point group $C_{s}^{d}$ in the third table of Fig. 10) were used, the two resonances would be cross-polarized. Tuning these three portions of the perturbation therefore allows independent tuning of each of the linewidths and the spacing of the two resonances in either a co-polarized or cross-polarized fashion. The co-polarization of the two previously degenerate resonances is unique to the $Hex_{K}$ lattice in the catalog, as the $E_{1}$ and $E_{2}$ modes of the hexagonal lattices (and $E$ modes in the square lattices) are only accessible in a cross-polarized fashion. Co-polarized split degenerate states are a unique feature of lattices with more than two atoms (such as the lattices containing four atoms seen in, for instance, Fig. B 1). In the present application, the freedom to have the pump wavelengths be co-polarized allows a single pulse (with bandwidth spanning the two resonances) as the pump.

To complete the demonstration of the advantages of the $Hex_{K}$ lattice for Terahertz generation, we compute a simple figure of merit related to the efficiency of this conversion (see, for instances, Refs. Lin et al. (2016); Lee et al. (2014)):

[TABLE]

where the bounds of integration are over the entire device and $\chi^{(3)}(x,y)$ is the spatially dependent third-order nonlinear susceptibility and the electric fields are normalized to the corresponding incident fields. Given the scale difference of $\lambda_{3}$ to a unit cell (i.e., $\lambda_{3}^{2}\gg A$), a reasonable approximation to this integral is that $E^{*}(\omega_{3})=F_{3}$ and $E^{*}(\omega_{4})=F_{4}$ are constants equal to the electric field enhancement due to the bowtie antenna. We may then integrate over a unit cell:

[TABLE]

where $F_{y}=E_{y}^{*}(\omega_{1})E_{y}(\omega_{2})G(x,y)$ and $G(x,y)=1$ where there is Silicon and is [math] where there is vacuum. That is, the figure of merit is proportional to the overlap integral of the two pumps within the Silicon portion of a unit cell.

The integrand may be calculated from the mode profiles taken from full-wave simulations of unit cell of the device. The $H_{z}$ component of the modes for the spectrum in Fig. 12(b) are shown in Figs. 12(f) and 12(g), corresponding to the choice in Fig. 12(c) of $R_{1}=0.38\mu m$. The integrand of Eq. (27) is shown in Fig. 12(h), as calculated from the $E_{y}$ components of the modes (shown in Figs. 12(i) and 12(j)) and the refractive index profile of the device. The numerical value for this case is $FoM/\chi_{Si}^{(3)}\approx 114\lvert F_{3}\rvert\lvert F_{4}\rvert$, meaning that for a modest enhancement of $\left\lvert F_{3}\right\rvert=\left\lvert F_{4}\right\rvert\approx 10$ by the bowtie antenna we will have a total enhancement of $FoM/\chi_{Si}^{(3)}\approx 10^{4}$ in efficiency. Figure 12 therefore demonstrates a platform to produce THz light from two infrared pumps taking advantage of large electric field enhancement at every frequency involved in the four-wave mixing process. The phase matching condition is guaranteed by the subwavelength scale of the device in the vertical direction, the resonance spacing is robustly controlled by the radius of the central pillar, and the resonant linewidths can be tuned largely independently by the successive degeneration from $C_{6v}$ to $C_{2v}$ to $C_{s}$. We note that there is some partial cancellation upon integration of the integrand of Eq. (27), but not complete cancellation. Future work could optimize the perturbations chosen such that this cancellation is minimized.

IV.3 Application: Mechanically Tunable Optical Lifetimes

Next, we remark on a type of periodic perturbation achieved by stretching or shearing a high symmetry lattice. Since the symmetry of the lattice is reduced, the symmetry-protected BICs may be excited. However, the condition on coupling (that is, Eq. (23)) still applies, and therefore this class of perturbation follows the same selection rules as the equivalent point group degeneration entries. For instance, by shearing an unperturbed $Sq_{\Gamma}$ lattice’s unit cell from a square into a rhombus, the space group is reduced to $p1$ and any mode at the $\Gamma$ point may now couple to free space with a strength related to the degree of shear. However, the polarization direction of the coupled planewave will be ill-defined in general, changing, for instance, with the degree of shear. (Recall that it is for this reason that the $p1$ entries are all specified as $any$, because no general comment can be made.)

Of more interest is stretching along a high symmetry axis, affording well-defined selection rules. This has been explored in Ref. Cui et al. (2012) for plasmonic heptamers arranged in a square lattice by degenerating the symmetry of the heptamer from $C_{6v}$ to $C_{2v}$ by stretching the substrate. Since the Fano resonance in the plasmonic heptamer is both (1) well-confined to a unit cell of the overall lattice, and (2) due to the coupling between plasmonic modes, analysis of the point group of the unit cell alone suffices to analyze the resonance. However, for a low loss, high Q-factor demonstration using dielectric structures, this analysis is insufficient because the coupling across unit cells of the array is integral to the presence of BICs. The catalog of selection rules derived above provides the necessary information for proper analysis in dielectric systems.

Inspection of the catalog of the square lattices (Fig. 9) reveals that square lattices afford no interesting cases: the only impact of a lattice deformation along a high symmetry axis is to split degeneracies, not introduce any new coupling. This is not true for the hexagonal lattices, however. Figure 13(a) shows a $Hex_{K}$ lattice with $p6mm$ space group (one pillar of the trimer has a larger radius than the others) on a stretchable substrate. In the unstretched case, the lattice has $C_{6v}$ symmetry, and the selection rules forbid coupling to any but the $E_{1}$ modes at normal incidence. However, inspection of the $cmm$ space group reveals that degeneration from $C_{6v}$ to $C_{2v}$ enables coupling to the $B_{1}$ and $B_{2}$ modes. Stretching the $Hex_{K}$ lattice with $p6mm$ space group along the $x$ axis also degenerates the point group from $C_{6v}$ to $C_{2v}$, and so ought to enable coupling to those modes to a degree controlled by the strength of the lattice deformation. Figure 13(b) depicts confirmation of this prediction via full-wave simulations near the $\text{TM}_{K,B_{1}}^{1,1}$ mode, showing redshift and a changing Q-factor as a function of deformation. Inspection of the $Hex_{K}$ catalog (or an analogous case in the $Hex_{\Gamma}$ lattice) therefore enables a low-loss dielectric-based flexible device platform with mechanically tunable resonant lifetime.

IV.4 Nonlocal Metasurfaces

The preceding applications demonstrated the utility of the catalog to guide device design using successive perturbations. In other words, the key to their design came from understanding how the final, lower order space group was constructed from the higher order space groups, which together form a parent-child relationship. In this section, we supplement the approach in the accompanying Letter Overvig et al. (2020) in which we use this principle to demonstrate how particular $p2$ space groups may be constructed from two parent space groups. We elucidate a geometric phase that is a consequence of the parent space groups exhibiting mode “twisting” of the sort shown in Fig. 11, and how if this geometric phase is spatially varied, we may realize devices with anomalous reflection and refraction only on resonance.

We begin by focusing on the relationship between three space groups in the $Sq_{M}$ lattice, shown in Fig. 14(a). Two parent space groups ($pmg$ and $pmm$) are shown on the left. These two space groups share no symmetries in common except two-fold rotations at the center of the Silicon pillars. Consequently, if the perturbations are added successively, the child space group will retain only these two-fold rotations, resulting in the $p2$ space group shown on the right in Fig. 14(a). This example of a parent-child relationship between higher order parents and lower order children is very general; the full hierarchy for the $Sq_{M}$ and $Sq_{X}$ lattice is reported in Fig. C 1. (Note: the hierarchy particular to the multi-wavelength metasurface introduced in the accompanying Letter Overvig et al. (2020) is reported in Fig. C 2)

Examining Fig.14(a) shows that the $pmm$ parent space group allows for coupling to a polarization angle $\phi=0\degree$ from the $x$ axis (that is, $x$ polarization), while the $pmg$ parent space group allows for coupling to a polarization angle $\phi=90\degree$ from the $x$ axis (that is, $y$ polarization). As seen in the inset of Fig.14(a), the child space group may be parameterized by an orientation angle, $\alpha$, that yields the $pmm$ parent space group when $\alpha=0\degree$, and the $pmg$ parent space group when $\alpha=45\degree$. In other words, as $\alpha$ varies continuously from $0\degree$ to $45\degree$, the corresponding polarization angle must vary from $0\degree$ to $90\degree$. The linear interpolation of this behavior is

[TABLE]

This form is highly reminiscent of the well-known geometric phase, $2\alpha$, which is introduced when the handedness of circularly polarization is flipped while light is scattered by an anisotropic scatterer oriented along the $\alpha$ direction. This similarity suggests studying a $p2$ space group under circularly polarized illumination, as shown in Fig. 14(b) for RCP light.

In particular, upon studying the phase of circularly polarized light exiting the device on both the reflection side and transmission side, we find that this system imparts a geometric phase to light exiting with the converted handedness (for RCP incidence, this is LCP in transmission, and RCP in reflection) that is twice the conventional geometric phase. As in conventional dichroic optical elements (e.g., a plasmonic bar antenna), two projections of the polarization basis are required to analyze the outgoing light, one from coupling into the element and the second from coupling out. We will consider each in turn.

First, only the component of free-space light that is linearly polarized in the $\phi$ direction completley couples to the mode. This light, which constitutes half of the power of the RCP incident light, is resonantly reflected, while orthogonally polarized light (at an angle $\theta=\phi+90\degree$) is transmitted. Decomposing the incident RCP light into two linearly polarized components, the $\phi$ component carries a phase $\Phi^{r}_{1}=\phi$, and the $\theta$ component carries a phase $\Phi^{t}_{1}=\phi+90\degree$. Since the $\phi$ direction is defined by the orientation angle $\alpha$ by Eq. 28, the resonantly reflected light is therefore associated with a phase $\Phi^{r}_{1}=2\alpha$ and the orthogonally polarized transmitted light has a phase $\Phi^{t}_{1}=2\alpha+90\degree$.

Second, the output light on resonance is linearly polarized, and can be decomposed into its constituent LCP and RCP components. These components have a geometric phase $\Phi^{r}_{2}=\mp\phi$ in reflection and $\Phi^{t}_{2}=\pm\theta$ in transmission, where the first sign corresponds to LCP, and the second corresponds to RCP. We can finally determine the total phases $\Phi=\Phi_{1}+\Phi_{2}$ of the LCP and RCP components in reflection and transmission:

[TABLE]

That is, RCP light in reflection and LCP light in transmission vary as $4\alpha$, while LCP light in reflection and RCP light in transmission are invariant to $\alpha$. The two factors of the geometric phase come from the two instances changing the basis for the polarization state: first from circular to linear (coupling into a single linear state), and second from linear to circular (decomposing into its constituent spins). In this case, the final value is twice the conventional geometric phase because the eigenpolarization (characterizing the projection bases) varies as $\phi=2\alpha$ compared to the conventional case of $\phi=\alpha$ (e.g., a plasmonic bar antenna oriented in-plane by an angle $\alpha$). We note that the form $\phi=2\alpha$ is not a general rule for quasi-BICs; for instance, the $p2$ space group in the $Hex_{K}$ lattice with a cross-motif follows $\phi=-4\alpha$ for the $B_{1}$ mode. This results in a geometric phase that is $-8\alpha$.

We next explore the physics and applications exploiting this geometric phase as a new degree of freedom. Because the geometric phase is completely controlled by $\alpha$, we may use the remaining geometric degrees of freedom of the unit cell to maintain a spatially constant resonant frequency across a device with a spatially varying geometric phase profile. In other words, as in conventional metasurface approaches, we may construct a library of geometries (“meta-units”) such that full phase coverage and constant amplitude are achieved. Then, by spatially arranging these meta-units, an output wavefront with a designer phase profile may be realized at the resonant frequency.

To confirm this approach, we construct such a meta-unit library targeting wavelengths in the telecommunications range. A meta-unit, seen in Fig. 14(b), is composed of a Silicon slab with two etched ellipses, which are identical but for a $90\degree$ rotation. The chosen thickness of the slab is $H=250nm$, and the lattice constant is $a=400nm$. The in-plane geometric parameters that are varied to construct the meta-unit library are the diameters along the semi-major axis, $D_{a}$, and semi-minor axis, $D_{b}$, and the orientation angle, $\alpha$. For simplicity, we keep $D_{b}$ constant, and vary $D_{a}$ and $\alpha$ so as to achieve full phase coverage with minimal shift in resonant wavelength.

A spectral map of reflectance, calculated by fullwave simulations, is shown in Fig. 14(c) illustrating a near constant resonant wavelength across the meta-unit library. The amplitude and phase of the reflected (transmitted) LCP and RCP components are recorded in Fig. 14(d) (Fig. 14(e)) at the operating wavelength, $\lambda_{op}=1.52\mu m$. The amplitudes of the LCP and RCP components are approximately equal (each representing roughly one quarter of the input power) and vary little across the meta-unit library. The small inequality is due to the presence of the substrate breaking the symmetry in the out-of-plane direction. The phase of the component with converted handedness (which, for reflection is RCP, and for transmission is LCP) varies across $2\pi$ as $\alpha$ varies across $90\degree$, and follows closely with the predicted $\Phi=4\alpha$ dependence (see Figs. 14(d,e)).

With the meta-unit library constructed, a wavefront with a spatially shaped phase profile within a narrow bandwidth near $\lambda_{op}$ may be realized. A common function is to linearly vary the output phase so as to create anomalous reflection and refraction. We may choose either to vary the phase profile in the same direction as the dimerization (the $x$ direction in Fig. 14(b)), or the orthogonal direction. We will begin with the former choice (see the accompanying Letter for the latter Overvig et al. (2020)).

Figure 15(a) shows a schematic of a device deflecting the component with the converted handedness at the resonant wavelength for RCP light normally incident from the substrate side. Figure 15(b) depicts the electric field on resonance calculated by fullwave simulations, overlaid with the geometry of the device (Fig. 15(e)). Figure 15(c,d) show the output polarization states, and phases of the LCP and RCP components at the reflection side and transmission side, respectively. In both cases, the output polarization is approximately linear across the device. As such, the phase of the signal with unconverted handedness is uniform while the phase of the signal with converted handedness varies across $4\pi$ as $\alpha$ varies over $180\degree$.

Figure 15(f) confirms that the resonance of the metasurface remains intact, despite the variance of geometry across the metasurface. However, a noticeable blueshift has occurred relative to the originally chosen $\lambda_{op}$. Nevertheless, at the resonant peak of the device, $\lambda_{dev}=1.46\mu m$, deflection to the second diffractive order occurs for signal with the converted handedness (Figs. 15(g-j)). A device with identical deflection angle is shown in the accompanying Letter Overvig et al. (2020), but with a phase gradient applied in the orthogonal direction to the dimerization direction. The blueshift is also present in that case, but significantly reduced.

The explanation for this blueshift comes from a unique feature of this metasurface: that the deflection of light with converted handedness is mediated by a supermode of the device. That is, unlike conventional metasurfaces, whose meta-units scatter light based on local resonances, this metasurface scatters light due to a global resonance (associated with a supermode) supported by many neighboring meta-units. To explore the physics here, we consider the dependence of a gradient resonant metasurface on the incident angle of the RCP light. Figure 16(a) schematically shows a device with a spatial phase gradient in a direction orthogonal to the dimerization direction, with light incident from the substrate at a set of angles $\theta_{x}$ (along the phase gradient) and $\theta_{y}$ (along the dimerization direction). Since the deflection only occurs on resonance, the resonant frequency follows some dispersion relation (i.e., the band structure). Figure 16(b) depicts the resonant frequency dispersion while varying $\theta_{y}$ from $-40\degree$ to $40\degree$ (corresponding, by Snell’s Law, to $\pm 68.7\degree$ in air). This mode is concave up, meaning that at higher in-plane momenta, a blueshift occurs.

To understand the blueshift at normal incidence, we must consider (1) the modes supported by the device and (2) the in-plane momentum of the resonant mode. First, these resonant modes exist in the device with a super-period of $P_{x}=8a$ in the phase gradient direction. We therefore consider all of the supermodes supported by a device with this super-period. Because the super-period is composed of perturbed versions of the same PCS, the supermodes will be well approximated by artificial Brillouin folding corresponding to period doubling a meta-unit three times. Figure 16(c) depicts such a process for the mode in question, showing the band of the unit cell (containing two ellipses) copied every integer multiple of the grating vector $k_{G}=2\pi/P_{x}$ in the $k_{x}=k_{0}sin(\theta_{x})$ direction. The supermodes present at normal incidence are the modes at $k_{x}=0$ in this diagram.

Second, upon coupling in, there is a spatially varying geometric phase, corresponding to twice the local rotation angle of the ellipses. The derivative of this spatial phase is equivalent to a $k$-vector,

[TABLE]

In other words, the resonant supermode is the resonant mode of a unperturbed lattice modulated by an in-plane wavevector in the $x$-direction. This corresponds to a supermode that is $k_{geo}$ away from the unperturbed $\Gamma$ point, $k_{x}=0$ (the flat part of the band). The supermode is highlighted red in Fig. 16(c). Consequently, as seen in Fig. 16(d), as $\theta_{x}$ is varied, the resonance follows the dispersion of the band as it existed in the unperturbed lattice, shifted by $-k_{geo}$. Notably, this also means that the resonant frequencies corresponding to the pair of incident momenta $k_{x}=0$ and $k_{x}=-2k_{geo}$ are identical. This is consistent with the requirements of reciprocity: these two momenta are the input and output momenta of the deflection process; reversing the output must yield the original input at all frequencies.

Finally, we comment on the achievable phase gradient limited by this angular dispersion. As encapsulated by Eq. 2, the component $k$-vectors involved with the resonance must be limited according to the linewidth and angular dispersion of the resonance in order to maintain large resonance visibility. For a resonant metasurface lens or hologram shaping an incident planewave, the range of output $k$-vectors must satisfy Eq. 2. In the cylindrical metasurface lens reported in the accompanying Letter Overvig et al. (2020), the resonance visibility is maintained despite the range of deflection angles in the $\theta_{x}$ direction across the device, which may be characterized by the numerical aperture, $NA$. We find that increasing the $NA$ gradually reduces the resonance visibility, but a substantial resonance visibility (a peak reflectance of $>40\%$) is still maintained at a high value of $NA=0.7$. This is not true for a cylindrical lens focusing in the $\theta_{y}$ direction, where $NA<0.1$ is required to retain appreciable resonance visibility, as the angular dispersion is large in the $\theta_{y}$ direction compared to that in the $\theta_{x}$ direction (comparing Figs. 16(b,d)). This is consistent with the constraint on Q-factor, band flatness, and spread in incident k-vector encapsulated by Eq. 2. In other words, it suggests that by including band structure engineering, a fully radially focusing resonant metasurface lens may be realized, and that we must generally take care to engineer the band structure of the unperturbed resonant metasurfaces before applying the perturbation, as laid out in the three-step process described in Sec. II.1.

Given the lack of impact on the non-resonant light waves, which may transmit with high efficiency irrespective of incident angle, we anticipate these resonantly deflecting and focusing metasurfaces to be of significant interest to augmented reality displays, which aim to superimpose a desired image on top of information transmitted through the glass originating from the external world. By further application of the principle of successive perturbations, we show in the accompanying Letter Overvig et al. (2020) that the single-wavelength resonant metasurfaces may be extended to multi-wavelength devices with independently tunable phase profiles. The hierarchy of the child space group constructed from eight parents is shown in Fig. C 2. The eight parents represent eight degrees of freedom to spatially and spectrally shape an incident wavefront: the Q-factors and polarization angles (i.e., geometric phases) of four modes with distinct symmetries may be controlled simultaneously. Notably, these eight degrees of freedom are in addition to the degrees of freedom present in the unperturbed lattice, which may be used to control the resonant frequencies and band curvatures of the desired modes. This degree of spatial and spectral control over an optical spectrum greatly surpasses the state-of-the-art, and is readily apparent from careful study of the catalog in conjunction with the design principle of successive perturbations.

V Summary

In summary, we derived the selection rules for Fano resonances due to quasi-bound states in the continuum supported by photonic crystal slabs. Targeting the high symmetry modes of both square and hexagonal lattices, we explored six lattices designed to bring each class of high symmetry mode into the continuum. We exhaustively reported the degenerated space groups due to in-plane perturbations that are compatible with these six lattices, and cataloged the selection rules in each case by applying principles of Group Theory to determine the free space polarization of the leaky portion of the perturbed modes.

Together with band structure engineering, the principles, approach, and results outlined here provide a high-level guide to designing compact photonic crystal slabs supporting sharp resonances: devices confining light in both space and time, and manufacturable by mature fabrication technologies. Future work will be well-guided by the rational design principles considered here to reduce the search space required to optimize a compact, resonant optical device. In particular, we showed that the band structure may be engineered in the unperturbed lattice before a periodic perturbation is applied to couple the targeted mode(s) to the desired free space polarization(s). We showed, here and in the accompanying Letter Overvig et al. (2020), that in addition to the degrees of freedom present in choosing the unperturbed lattice, a series of successive perturbations may realize multifunctional control of the resonances (up to eight parameters at once). The insights of the catalog of selection rules produced by Group Theory arguments have straightforwardly motivated novel devices, such as polarization independent planar optical modulators, Terahertz generation in photonic crystal slabs with lifted degeneracies, devices with mechanically tunable optical lifetimes, and a novel class of metasurfaces that uses two factors of a geometric phase to spatially shape a resonant wavefront. We therefore believe that careful understanding and examination of the patterns and features of the selection rules represent a fruitful launching point for future efforts.

Acknowledgements.

The authors thank J. Lee-Thorpe and M. Weinstein, and Andrea Alù for helpful discussions. The work was supported by the Defense Advanced Research Projects Agency (grant nos. D15AP00111 and HR0011-17-2-0017), the National Science Foundation (grant nos. ECCS-2004685 and QII-TAQS-1936359), and the Air Force Office of Scientific Research (grant no. FA9550-14-1-0389). A.C.O. acknowledges support from the NSF IGERT program (grant no. DGE-1069240). S.C.M acknowledges support from the NSF Graduate Research Fellowship Program (grant no. DGE-1644869). M.J.C. acknowledges support from USD/R&E (The Under Secretary of Defense-Research and Engineering), National Defense Education Program (NDEP)/BA-1, Basic Research.

Appendix A Group Theory Tables

For ease of reference, the character tables of all relevant point groups are reported in Fig.A 1(a). The left column of tables in Fig. A 1(a) contains the point groups compatible with the square and rectangle lattices, and the right contains those compatible with hexagonal lattices. Figure A 1(b) summarizes the subgroups of each of the point groups shown in Fig. A 1(a). This prescribes the necessary components of the symmetry degeneration tables, shown in Fig. A 1(c), which track how higher symmetry modes (irreducible representations) degenerate into lower groups. That is, reference to the symmetry degeneration tables provides the answers to how a higher symmetry mode would be named in a lower order symmetry group (for instance, $B_{2}$ in $C_{4v}$ would be called $A_{2}$ in $C_{2v}^{v}$).

Next, Fig. A 2 provides the Group Theory tables helpful for determining the selection rules through the direct product approach. Figure A 2(a) provides the direct product table for the $C_{6v}$ point group, and Fig. A 2(b) provides the same for the $C_{4v}$ point group. The direct product tables for the lower order point groups are a subset of these. For instance, for $C_{2v}$, Fig. A 2(b) may be used excluding the final row and column. Lastly, Fig. A 2(c) contains the irreducible representations for the partial derivative operators relevant for direct products such as in Eq. (23). This tracks how free space polarizations transform in each lattice type.

Appendix B Additional Lattices

The six lattices catalogd in Figs. 9 and 10 were chosen because they access the six high symmetry modes in the simplest way. For instance, the $Sq_{M}$ lattice accesses only the $\Gamma$ and $M$ (but not the $X$) modes of a square lattice, while the $Sq_{X}$ accesses the $\Gamma$ and $X$ (but not the $M$) modes. However, if desired, it is possible to access all three modes in a single lattice by period doubling in both lattice directions. Figure B 1 shows this lattice, called $Sq_{MX}$, depicting the real space and FBZ. This lattice is a ”quadromer”, having four atoms per unit cell, and therefore has four times the modes compared to the unperturbed case. In particular, it has the $\Gamma$ and $X$ modes, as well as two copies of the $M$ modes (one from each $M_{x}$ and $M_{y}$), which mix at the $\Gamma$ point in a similar way to the $K$ modes of the $Hex_{K}$ lattice.

Also pictured in Fig. B 1 are the $Hex_{MMM}$ and $Hex_{M\Lambda}$, both of which are examples of quadromer lattices. $Hex_{MMM}$ is still of the hexagonal lattice family, while $Hex_{M\Lambda}$ (much like $Hex_{M}$) is rectangular. $Hex_{MMM}$ contains three copies of the $M$ modes, which will mix at the $\Gamma$ point. $Hex_{M\Lambda}$, on the other hand, accesses a unique set of modes at the $\Lambda$ point in the unperturbed FBZ (see the last panel in Fig. B 1), which have the point group $C_{s}^{d}$ (that is, they are either symmetric or anti-symmetric about the $x$ axis).

The lattices shown in Fig B 1 are by no means the only additional lattices that may be explored. Instead, they serve as an example of the next few lattices in the infinite list of lattices ordered by number of atoms in the perturbed unit cell. The lattices in this list are generally increasingly complicated, but the same approach outlined in Sec. III.3 may be applied to determined the selection rules if desired.

Appendix C Space Group Hierarchy

A feature of the $Sq_{X}$ lattice is that there are two points in the unperturbed lattice with $C_{4v}$ about which to apply the dimerizing perturbation. This has the consequence that there are two $p4g$ groups in the catalog of the $Sq_{X}$ lattice, each with a different high symmetry point in common with the unperturbed lattice (likewise for the two $p4m$ groups). Similarly, there are two equivalent points with at least $C_{2v}$ symmetry that may be chosen while perturbing to a $Sq_{M}$ lattice (and equivalently, the $Hex_{M}$ lattice), producing an analogous set of paired perturbations. Note that the $Sq_{\Gamma}$ lattice (and equivalently, the $Hex_{\Gamma}$ lattice) has no such feature. The $Hex_{K}$ lattice has no such pairing of perturbations because there is only one point in the unperturbed lattice with $C_{6v}$ symmetry.

An interesting feature of the paired sets of perturbation is the relationship of the symmetries those members have in common with each other. By adding the two perturbations together, a space group of lower symmetry is created with only those symmetries the two parent space groups had in common. In the words of crystallography, this process is finding a translationengleiche subgroup in common (and finding the $C_{4v}$ parent space groups is finding klassengleiche subgroups of the unperturbed lattice). This resulting child space group must, naturally, be a member of the catalog. Such a relationship may be studied for all of the space groups in the lattices with this pairing of perturbations. The resulting relationships are seen in Figure C 1, showing complete sets of relationships between higher and lower space groups. Note that for the $Sq_{X}$ lattice, one of the $cmm$ space groups has its highest symmetry point in common with the $C_{2v}$ symmetry point of the unperturbed lattice. It may therefore not be made by adding higher order space groups. Likewise, it is well known to crystallographers that the $pmg$ space group has no parent (translationengleiche supergroup), and therefore also stands alone. However, both of these “parent-less” space groups share symmetries with the other members of the $Sq_{X}$ lattice, and therefore may be combined to create space groups of lower symmetry.

We note that while these relationships do not prove that the catalog is exhaustive, the closed, consistent system is highly suggestive that it is. Indeed, in an earlier version of this manuscript, one of the $p4g$ space groups of the $Sq_{X}$ catalog was omitted, and its existence and selection rules were predicted while attempting constructing Figure C 1(b): a fourth space group with $C_{4v}$ symmetry was needed to produce a consistent hierarchy.

Finally, for completeness we report in Figure C 2 the entire hierarchy enabling control of the Q-factor and polarization angle for four modes simultaneously and independently, used in the accompanying Letter Overvig et al. (2020). The resulting lattice is $Sq_{MX}$ (refer to Figure B 1) constructed from two $Sq_{M}$ lattices, one dimerized in the $x$ direction, whose modes are referred to as $M_{x}$, and the other in the $y$ direction, whose modes are referred to as $M_{y}$.