Topological insulators (TIs) were originally discovered in condensed-matter systems [1,2], and then further explored in various disciplines such as electronics [3], mechanics [4], and phononics [5]. In particular, we have recently witnessed the rise of topological photonics, which endows light with intriguing topological properties [6–8]. For example, topological edge states, the boundary states existing on the surfaces of TIs, are introduced as guide modes in photonic circuits, since they can support robust light-wave propagation immune to backscattering and perturbations, having potential applications in elevating transmission quality and efficiency of optical signals [9,10]. More recently, a new member of TIs, higher-order TIs (HOTIs), has drawn increasing attention due to the lower-dimensional gapless boundary states. For instance, two-dimensional (2D) second-order TIs (SOTIs) have a 1D edge and a 0D corner states. Distinct from the propagating nature of edge states, corner states are highly localized in the corner and have been identified in phononic, photonic, and circuit systems [11–13]. Such a feature also provides a novel approach to trap light in photonic crystal (PhC) platforms [14–16], and has been applied for exploring cavity quantum electrodynamics [17] as well as high-performance lasing [18,19].

At the same time, strong coupling between light and matter has attracted much interest [20,21]. When the coherent energy exchange rate between photons and matter excitations exceeds their initial decay rates, they come to the strong-coupling regime and form the half-light, half-matter quasi-particles known as polaritons [22,23]. Up to now, strong coupling accompanied with polaritons has opened up new avenues for not only light manipulation but also modification of material properties and chemical reactions [24–26]. Especially, when photons gather in topological states and strongly couple to matter, they give rise to topological polaritons [27,28]. For example, topological exciton-polaritons have emerged as a promising platform linking topological photonics and matter. Up to now, topological exciton-polaritons have been experimentally realized in GaAs microcavities and PhC slabs assembled with 2D materials, paving the way towards topological polaritonic devices [29–32].

Another remarkable research topic in topological photonics is endowing the topological states with the capability of on-demand dynamically reconfiguration, since it is essential in practical applications [33]. Nowadays, various schemes have been employed to realize tunable topological states, including mechanical control [34,35], liquid crystal reorientation [36,37], thermal-induced state transition of phase-change material [38], optical control based on free-carrier excitation in silicon [39], and modifying the defect strength [40]. Especially, as a typical 2D material, graphene also has potential applications in electrically reconfigurable topological devices, since its Fermi level, which is related to its material properties, can be readily tuned by externally applied voltage [41,42].

Recent advances have revealed that topological states can exist as leaky states, namely, couple to the radiative continuum, in metasurfaces. Such non-Hermitian systems can be excited by far-field light sources and serve as a versatile platform for various applications such as topological invariant measurement [43], nonlinear effect enhancement [44], and topological polariton generation [31,32]. Thus, it is natural to wonder if it is feasible to introduce localized topological states to this system, since metasurfaces’ building blocks (meta-atoms) usually support localized modes [45,46]. In addition, on the side of strong coupling, there is growing interest in exploring the interaction between photons and lattice vibrations in polar materials [47–49]. The accompanied phonon polaritons, with low-loss and energy-concentration features, may play an essential role in terahertz and mid-infrared polaritonic devices in the future [50,51]. More recently, one work has implemented the robust and unidirectional transport of the topological phonon polaritons by the strong coupling between the helical edge state and phonons in the mid-infrared [52]. However, the topological phonon polaritons mediated by corner states, which could provide a robust approach to localize phonons, have not been explored yet. Furthermore, from the application’s perspective, active control towards topological states and strong coupling is also indispensable.

In this work, we demonstrate a supercell metasurface that supports localized corner states and investigate its interaction with monolayer hexagonal boron nitride (hBN) and graphene. The metasurface comprises periodically arranged PhC-slab SOTIs. Consequently, the corner states in supercells go beyond the light line due to the shrinking of the first Brillouin zone (FBZ) and can be stimulated by a far-field source. Via the plane-wave expansion (PWE) method, we show that the topological transition still exists even in the absence of the bandgap in unit cells, and nondegenerate multipolar corner states emerge due to the inter- and intrasupercell coupling effects. The full-wave simulation reveals that the two dipolar corner states, as well as their superposition state, can be collectively excited by a far-field plane-wave source in the mid-infrared regime, showing a sharp dip in the transmission spectrum. By integrating monolayer hBN with the metasurface, we further investigate the strong photon–phonon coupling between corner state and phonon in hBN via spectral dynamics, where the Rabi splitting and anticrossing features can confirm the occurrence of the strong coupling. The influences from material properties, as well as positions of hBN, are studied in detail. Moreover, the defect is introduced to test the topological protection of this structure. Finally, we show that by assembling monolayer graphene onto the metasurface, the resonant frequency of the corner state can be sufficiently tuned via modifying the Fermi level of the graphene so as to realize an electrically reconfigurable corner state and strong coupling. This work introduces a unique platform to the realm of topological photonics, and the resulting features could be applied to large-area topological phonon polaritonic devices and reconfigurable topological devices in the mid-infrared region.

The geometrical configuration of a supercell metasurface is shown in Fig. 1. It is composed of periodically arranged PhC-slab supercells along the $x$ and $y$ directions, each of which consists of a finite number of square unit cells, as defined in Fig. 1. The lattice constants of the unit cell and supercell are $a$ and $5a$, respectively. Here we set $a=2.03\mathrm{\mu m}$ so that the corner states will be located in the mid-infrared region. The height of the metasurface is $H=4.245\mathrm{\mu m}$. The material of the metasurface is chosen as silicon, whose dispersionless dielectric constant is $\epsilon =11.7$ in the mid-infrared region [53]. As presented in Fig. 1, there is a square region filled with nontrivial unit cells (highlighted in blue) surrounded by trivial unit cells (highlighted in gray) in each supercell, where $m$ represents the number of trivial rows between neighboring nontrivial regions and $n$ denotes the number of nontrivial rows in one supercell. Here, the “nontrivial” and “trivial” refer to the topological nontrivial and trivial phases that unit cells possess, respectively, which will be discussed in the following section.

To reveal the topological properties of different 3D unit cells, first, we employ the PWE method to compute photonic band structures of PhC-slabs with these unit cells; the calculated results are summarized in Fig. 2. Notably, here only the transverse electric (TE) modes are considered in the calculation, because the bandgaps, which play an important role in the topological transition, are favored for TE modes in such connected periodic structures [54].

The primary unit cell, as shown in Fig. 2(a), contains four identical square air holes with side length $l=0.34a$, clustering together at the center of the unit cell. Our model emulates the 2D Su–Schrieffer–Heeger (SSH) model, whose Hamiltonian is as follows [14]: $$H=\left(\begin{array}{cccc}0& h_{12}& h_{13}& 0\\ h_{12}^{*}& 0& 0& h_{24}\\ h_{13}^{*}& 0& 0& h_{34}\\ 0& h_{24}^{*}& h_{34}^{*}& 0\end{array}\right),$$(1)where $h_{12}=h_{34}=t_a+t_b\exp (i{k}_x)$, and $h_{13}=h_{24}=t_a+t_b\exp (-i{k}_y)$. Here, $\mathit{k}=(k_x,k_y)$ is the wave vector, $t_a$ and $t_b$ are intra- and intercell hopping among atoms in the SSH model. In our PhC unit cell, the four holes play the role of the atoms in the SSH model, and shifting the holes along the diagonal of the unit cell will tune the intra- and intercell optical coupling, serving as the analog of hopping modulation in the SSH model.

Following the scheme mentioned above, we move the four holes outwardly along the diagonal of the unit cell, as shown in Figs. 2(a)–2(c), and observe the evolution of the corresponding band structures in FBZ, as presented in Figs. 2(d)–2(f) (only bands below the light line are included). In Figs. 2(d) and 2(f), there is no bandgap. However, we find two bands, highlighted in blue and red, respectively, that resemble the first two bands of the 2D counterpart of the unit cells, i.e., the 2D unit cells with infinite height. The $H_z$ field distributions of these two bands at symmetry points are similar to those of 2D unit cells (see Appendix A). The only difference is that the bandgap, which exists in band structures of 2D unit cells, is populated with several emerged bands (highlighted in gray) in the band structures of 3D unit cells. Furthermore, we find that the two bands are fundamental state bands with no nodal plane along the $z$ axis, whereas those emerged bands are higher-order state bands, which possess more than one nodal plane along the $z$ axis (see Appendix B). Previous studies have revealed the existence of these higher-order state bands, which are known as photonic subbands [55,56]. These states result from the confinement of optical fields within the slab, and can be regarded as the higher-order guide modes in a uniform slab distorted by in-plane hole lattices [56]. When the PhC-slab is thin, these states are located at higher frequencies than fundamental states due to the lower effective refractive index ($n_{\mathrm{eff}}$). Lifting the height of the slab will increase the $n_{\mathrm{eff}}$ of all states and redshift them. However, the $n_{\mathrm{eff}}$ of higher-order states are more sensitive to such change, so they will gradually fall into and cancel the fundamental bandgap [56–58], leading to the consequences shown in Figs. 2(d) and 2(f).

Next, we focus on the two fundamental state bands of 3D unit cells. The field distributions of $H_z$ of the two bands at the $\mathit{X}$ point are shown in the insets of Figs. 2(d) and 2(f), obtained in the $x-y$ plane crossing through the middle of the unit cell and perpendicular to the $z$ axis. Note that $H_z$ field profiles in the same $x-y$ cross section will be shown in the following section unless stated otherwise. In Fig. 2(d), the $H_z$ of the first band has even parity (represented by the symbol + in the inset, and indicated in blue band color), whereas another band has odd parity (represented by the symbol – in the inset, and indicated in red band color). As the four holes are moved outwardly along the diagonal until they occupy the four quadrants of the unit cell symmetrically, as shown in Fig. 2(b), the two bands are degenerate along the $X-M$ crystal symmetry axis, as depicted in Fig. 2(e). When the four holes continue moving outwards to the edges of the unit cell as shown in Fig. 2(c), the corresponding band structures are depicted in Fig. 2(f), which are identical to those in Fig. 2(d). Nevertheless, the parity of $H_z$ in Fig. 2(f) displays the opposite behavior. Such parity inversion reveals that even when the emerged higher-order state bands populate the bandgap, they will not interfere with the behavior of the two fundamental state bands. The two bands still have been inverted, accompanied by the topological phase transition. This topological phase transition is characterized by the 2D Zak phase, which is defined as follows [14]: $${\theta }_{\alpha }=2\pi P_{\alpha }=\frac{1}{2\pi }{\int }_{\mathrm{FBZ}}\mathrm{Tr}[A_{\alpha }(\mathit{k})]\mathrm{d}\mathit{k},\alpha =x,y,$$(2)Here, $\mathit{P}=(P_x,P_y)$ represents the 2D polarization, and $A_{\alpha }(\mathit{k})=i⟨{\psi }_{\mathit{k}}|{\partial }_{k_{\alpha }}|{\psi }_{\mathit{k}}⟩$ is the Berry connection. In the 2D unit cell, ${|{\psi }_{\mathit{k}}⟩}_{2\mathrm{D}}$ is the eigenstates of all the bands below the bandgap, namely, the eigenstates of the first band. Since the band inversion of the 3D unit cell is very similar to that in the 2D case, here, we define $|{\psi }_{\mathit{k}}⟩$ as the eigenstates of the first fundamental state band. The calculated results show that the unit cell in Fig. 2(a) has a 2D Zak phase (0, 0), meaning that it possesses a trivial topological phase. In contrast, the Zak phase is ($\pi$, $\pi$) for the unit cell in Fig. 2(c), indicating that it has a nontrivial topological phase. The details of this calculation are given in Appendix C.

Previous studies of SOTIs have realized topological corner states in the finite square PhC-slabs, comprising nontrivial unit cells surrounded by trivial unit cells [14,15,19]. There are four multipolar corner states in the bandgap with degenerate eigenfrequencies, namely, monopole, dipole I, dipole II, and quadrupole states [19]. Meanwhile, these corner states are protected by nontrivial topology, and characterized by the topological corner charge $Q_c=4P_x{P}_y=1$, where $P_x$ and $P_y$ are 2D polarization of the nontrivial unit cell [14]. Differently, in this work, we extend SOTIs to be periodic along the $x$ and $y$ axes, so that each SOTI acts as a supercell and a metasurface is formed, as shown in Fig. 1. When the number of trivial and nontrivial rows is $m=2$ and $n=3$, respectively, we obtain four corner states at $\mathrm{\Gamma }$ point of the metasurface, of which eigenfrequencies are in the gap between the first two fundamental state bands of the unit cells. The corresponding $H_z$ profiles are presented in the insets of Fig. 3(a). It is shown that the fields of these coupled corner states are highly localized at the corners of each supercell, and the distributions are consistent with the definition of multipolar corner states in Ref. [19]. At $\mathrm{\Gamma }$ point, the dipole I and dipole II states have degenerate eigenfrequencies. However, monopole and quadrupole states are nondegenerate, a feature that is distinct from the conventional corner states in an isolated SOTI. Such periodically coupled topological states have also been reported in previous works [43,59].

Such nondegeneracy mainly results from the interaction of the optical fields at corners, namely, the inter- (between neighboring supercells) and intrasupercell (within the same supercell) optical coupling among the corner states [19,59]. To reveal this point, first, we tune $m$ from 1 to 5 and keep $n$ at 5, and thus the intersupercell coupling strength is gradually weakened while the intrasupercell coupling strength maintains unchanged. As $m$ increases, the eigenfrequencies of multipolar corner states at $\mathrm{\Gamma }$ point converge asymptotically to a common value that corresponds to the eigenfrequency of the degenerate corner states in an isolated SOTI, as shown in Fig. 3(b). The same convergence manner can be found if $n$ is varied from 2 to 6 while $m$ is fixed at 4, as presented in Fig. 3(c). These phenomena indicate that these corner states will become degenerate if the coupling effects vanish, proving that the coupling effects will break the degeneracy of the corner states.

In particular, the periodicity of the metasurface endows the corner states with band structures, as shown in Fig. 3(a), which are given along the $k_x$ direction. Due to a lack of bandgap, these bands are embedded in the surrounding higher-order bulk bands that are not shown here, but do not hybridize with them. It is shown that these bands are symmetric to $\mathrm{\Gamma }$ point, and display strong frequency dispersion. In addition, the sign of the group velocity for the monopole and dipole II bands is opposite to that for dipole I and quadrupole bands, which is explained by the opposite convexity of the corresponding bands. Interestingly, the degeneracy of the two dipolar corner states can be broken upon the variation of $k_x$.

In addition, the corner states in this metasurface are far beyond the light line because the metasurface, whose lattice constant is $5a$, possesses shrinking FBZ compared with the system with lattice constant $a$, which will lead to the shrinking of the light cone. Thus, these corner states come to the continuum, become radiative, and couple to the external excitations. To reveal this, we perform a full-wave simulation based on the finite-difference time-domain (FDTD) method. In the simulated domain, the supercell is modeled with periodic conditions in the $\pm x$ and $\pm y$ directions and perfectly matched layers (PMLs) in the $\pm z$ directions. The far-field excitation is introduced by an $x$-polarized plane wave source at the top of the supercell, propagating along $-z$ direction, as shown in Fig. 1. One 2D power monitor across the $x-y$ plane is placed under the metasurface, and the total transmitted power $P_t$ is collected by integrating the power across the surface of the monitor. The normalized transmission spectrum is calculated by $T(\omega )=P_t(\omega )/P_s(\omega )$, where $P_s$ represents the power of the light source.

The transmission spectra of the metasurface with different height $H$ under the normal incidence are shown in Fig. 3(d). When $H=4.245\mathrm{\mu m}$, there is a transmission dip at $1387.21{\mathrm{cm}}^{-1}$ ($0.282a/\lambda$), which corresponds to the eigenfrequencies of dipole I and dipole II states in Fig. 3(a) at $\mathrm{\Gamma }$ point. It is a leaky state (resonant state) that is efficiently excited by the external stimulation, and the interference between the resonant state and the incoming plane wave gives rise to the Fano profile of the transmission spectrum [60], $$T_{\mathrm{Fano}}={\left|a_1+i{a}_2+\frac{b}{\omega -{\omega }_0+i\gamma }\right|}^2,$$(3)where $a_1$, $a_2$, and $b$ represent real constant parameters, ${\omega }_0$ represents the resonant frequency, and $\gamma$ represents the damping rate of the resonance. By fitting the transmission spectrum with this Fano formula, the quality ($Q$) factor can be calculated following the relation $Q={\omega }_0/2\gamma$. Here, the $Q$ factor of this resonance is 3022. In addition, the resonant frequency of the state can be slightly tuned by changing the height of the metasurface, as indicated in Fig. 3(d).

We further confirm that this resonant state is a dipolar corner state by investigating the field profile at the resonant frequency, as presented in the insets of Fig. 3(d). The $H_z$ distribution in the $x-y$ cross section shows that the field is localized at four corners of the supercell, which is the superposition state of dipole I and dipole II states. The reason for the emergence of this state is that the eigenfrequencies of dipole I and dipole II states degenerate at $\mathrm{\Gamma }$ point, and the field profiles of the two states can match the spatial symmetry of the incident plane wave polarization simultaneously, which is odd under 180° rotation around the $z$ axis. Especially, by rotating the polarization direction of the plane wave, we also successfully excite the two dipolar states separately (see Appendix D). In addition, the $H_z$ distribution in the $x-z$ cross section of the dipolar corner state shows that the field is primarily confined in the metasurface. The other two corner states, monopole and quadrupole states, however, cannot be excited owing to the symmetry mismatch, namely, zero overlap between their field profiles and incident waves. Since these two states are localized states beyond the light line, and do not couple to the continuum, they exist as bound states in the continuum (BICs).

What is more, the collective nature of the corner states in metasurface is revealed by investigating the responses of finite SOTI square arrays with different sizes of $N\times N$, where $N^2$ represents the total number of supercells. The normalized scattering cross sections and extracted $Q$ factors of these finite arrays are plotted in Figs. 3(e) and 3(f), respectively. The $Q$ factor of the dipolar corner state in an isolated supercell is only 152, as marked in Fig. 3(f). However, as the size of the finite array increases, the generated scattering cross section becomes narrower, and the peak approaches the resonant frequency of the infinite array, as shown in Fig. 3(e). Meanwhile, the corresponding $Q$ factor rises exponentially toward the $Q$ factor of the infinite array, as illustrated in Fig. 3(f). These features indicate a collective response of the metasurface; in other words, the response is not a sum of the individual supercell responses [61–63]. Such collective behavior arises from the coherent coupling between neighboring supercells, which can lead to the synchronous oscillations of the corner states in all supercells. The size-dependent $Q$ factor results from the truncation of the collective behavior caused by the finiteness of the array and will reach the level of the infinite array when the array is sufficiently large. Previous studies have revealed that such collective states can emit light with high spatial coherence, which can be applied for large-area lasing [64].

In this section, we present the strong coupling effect between dipolar corner states in an infinite array and lattice vibrations. The lattice vibrations, namely, phonons, are introduced by covering a monolayer hBN film (called hBN for simplicity in the following) onto the metasurface, as presented in the inset of Fig. 4(a). The surface conductivity model of hBN is as follows [65]: $${\sigma }_{\mathrm{hBN}}(\omega )=\frac{-4i{\epsilon }_0\omega {\omega }_{\mathrm{TO}}{v}_g}{{\omega }_{\mathrm{TO}}^2-{\omega }^2-i\omega {\gamma }_p}.$$(4)Here, ${\epsilon }_0$ is the vacuum permittivity, and $v_g$ is the phonon group velocity with the value of $1.2\times {10}^{-4}c$, where $c$ is the light speed in vacuum. Due to the absence of splitting between longitudinal and transverse optical (LO and TO) phonons, hBN possesses degenerate optical phonon frequency ${\omega }_{\mathrm{TO}}$ at $1387{\mathrm{cm}}^{-1}$ [66]. The intrinsic damping rate of phonon, represented by ${\gamma }_p$, can be artificially regulated by the isotopic enrichment of hBN [67]. Here, ${\gamma }_p$ is primarily set as $2{\mathrm{cm}}^{-1}$; different values will be discussed in this section.

The transmission spectrum of the suspended hBN under the plane-wave excitation is shown by the blue curve in Fig. 4(a), and there is a transmission dip at $1387{\mathrm{cm}}^{-1}$, which indicates the existence of the phonon. For comparison, the transmission spectrum of the dipolar corner state in the metasurface with $H=4.245\mathrm{\mu m}$ is also plotted in Fig. 4(a) as the red curve. The spectral overlap of the corner state and the phonon can facilitate coherent coupling between them.

With hBN laid on the top of the metasurface, a hybrid silicon–hBN metasurface is constructed, and a spectral splitting, namely, Rabi splitting, is observed in the transmission spectrum, as shown by the yellow curve in Fig. 4(a). There are two emerged transmission dips at $1385.09{\mathrm{cm}}^{-1}$ and $1389.02{\mathrm{cm}}^{-1}$, which indicates the formation of two new states. As depicted in Fig. 4(b), the two new states are lower-energy (${\omega }_{-}$) and upper-energy (${\omega }_{+}$) hybridized photon–phonon states, namely, topological phonon polariton states, which are the result of strong light–matter interaction between the photonic corner state and the phonon in hBN. A transmission dip at $1387.06{\mathrm{cm}}^{-1}$ still exists in the spectrum, representing nonhybridized photons and phonons remaining in the system [68]. In addition, the $|\mathit{H}|$ distributions of the two polariton states still maintain the profiles of the corner state, as presented in Fig. 4(c). The condition of strong coupling is as follows [20]: $$C=\frac{{\mathrm{\Omega }}^2}{\frac{{\mathrm{\Gamma }}_c^2}{2}+\frac{{\mathrm{\Gamma }}_p^2}{2}}>1,$$(5)where $C$ is the strong-coupling factor, $\mathrm{\Omega }$ is the Rabi splitting, and ${\mathrm{\Gamma }}_c$ and ${\mathrm{\Gamma }}_p$ are linewidths of the corner state and phonon, respectively, indicating the total damping rate induced by both radiative and intrinsic losses, and can be extracted from the corresponding transmission spectrum by fitting the Fano formula. Here, we get $\mathrm{\Omega }=3.9{\mathrm{cm}}^{-1}$, ${\mathrm{\Gamma }}_c=0.46{\mathrm{cm}}^{-1}$, and ${\mathrm{\Gamma }}_p=2.35{\mathrm{cm}}^{-1}$; thus, $C=5.31>1$, and the strong-coupling condition is fully fulfilled.

To analytically describe the coupling of the corner state and phonon, the coupled harmonic oscillators model (CHOM) is employed, in which the corner state and phonon are regarded as two coupled harmonic oscillators [69–71]. The two eigenfrequencies ${\omega }_{\pm }$ in this model, representing the two new polariton states, can be described as follows: $${\omega }_{\pm }+\frac{i{\mathrm{\Gamma }}_{\pm }}{2}=\frac{{\omega }_c+{\omega }_p}{2}-i\frac{{\mathrm{\Gamma }}_c+{\mathrm{\Gamma }}_p}{4}\pm \frac{1}{2}\sqrt{4g^2+{\left({\omega }_c-{\omega }_p-i\frac{{\mathrm{\Gamma }}_c-{\mathrm{\Gamma }}_p}{2}\right)}^2},$$(6)where ${\omega }_c$ (${\omega }_p$) represents the resonant frequency of the corner state (phonon), and $g$ is the coupling strength equal to $\mathrm{\Omega }/2$. Furthermore, we slightly tune the resonant frequency of the corner state ${\omega }_c$ by changing the height $H$; the corresponding transmission spectra are depicted in Fig. 4(d). An anticrossing behavior, another typical feature of strong coupling, is observed between the two polariton branches in Fig. 4(d), further proving that it is indeed in the strong-coupling regime. The ${\omega }_{\pm }$ as a function of ${\omega }_c$ are plotted in Fig. 4(e), showing that the numerical results can be fitted well by the analytical model where $g$ is treated as a free-fitting parameter. Additionally, as ${\omega }_c$ shifts, the weight of the corner state and phonon in upper and lower branches can be calculated as follows [72]: $${|\alpha |}^2=\frac{1}{2}(1\pm \frac{\mathrm{\Delta }}{\sqrt{{\mathrm{\Delta }}^2+4g^2}}),{|\beta |}^2=\frac{1}{2}(1\mp \frac{\mathrm{\Delta }}{\sqrt{{\mathrm{\Delta }}^2+4g^2}}),$$(7)where $\mathrm{\Delta }$ is the frequency detuning of ${\omega }_c$ relative to the primitive value $1387.21{\mathrm{cm}}^{-1}$; the corresponding results are shown in Figs. 4(f) and 4(g), respectively. It is shown that in the upper branch, the fraction of the corner state increases and the fraction of the phonon decreases as the detuning varies from positive to negative, as shown in Fig. 4(f), whereas the lower branch possesses the opposite behavior, as shown in Fig. 4(g). When $\mathrm{\Delta }=0$, the corner state and phonon contribute equally to the two polariton states.

To gain deeper physical insights into the strong coupling, the hBN is moved from the top ($z=H/2$) to the middle of the metasurface ($z=0$). The corresponding transmission spectra under different $H$ are depicted in Fig. 5(a), and the ${\omega }_{\pm }$ as a function of ${\omega }_c$ are plotted in Fig. 5(b). Interestingly, a Rabi splitting $\mathrm{\Omega }=7.24{\mathrm{cm}}^{-1}$ is achieved, and the anticrossing feature becomes more significant, as indicated in Figs. 5(a) and 5(b). Evidently, the strong coupling is further enhanced when the hBN comes to the middle of the metasurface. This is because the field of the corner state is concentrated more in the middle, as indicated in the inset of Fig. 3(d), and strong coupling will benefit from such local field enhancement [21]. In addition, the spatial overlap between the field and hBN increases, and more phonons will participate in the interaction, which will also contribute to the enhancement of the strong coupling.

The material properties of the hBN also play an important role in this strong coupling. To study this, we tune the phonon damping rate ${\gamma }_p$ of the hBN from 1 to $5{\mathrm{cm}}^{-1}$, meanwhile adjusting the position of the hBN from the top to the middle of the metasurface. The induced Rabi splitting $\mathrm{\Omega }$ and strong-coupling factor $C$ are plotted in Figs. 5(c) and 5(d), respectively. Notably, the missing data points represent the vanishing of Rabi splitting. As shown in Fig. 5(c), the higher phonon damping rate will lead to the decline of the Rabi splitting, provided hBN is at the same position. It is because a higher phonon damping rate means a shorter lifetime of the phonon, which will reduce the interaction time between photons and phonons. It is also shown that moving the hBN from the top to the middle will gradually enhance the Rabi splitting due to the field enhancement, and Rabi splitting is always visible for hBN at any position when ${\gamma }_p$ is within $1–3{\mathrm{cm}}^{-1}$. Notably, Fig. 5(d) indicates that the strong-coupling condition is fulfilled for any case plotted in Fig. 5(c), showing that such a hybrid metasurface can serve as an ideal platform for strong coupling. In addition, a similar evolution of $\mathrm{\Omega }$ and $C$ can be expected for the position of hBN from $z=0$ to $z=-H/2$ due to the symmetry of the field distribution along the $z$ axis.

Moreover, hBN may be etched together with the metasurface in practical fabrication. Thus, we further investigate the strong coupling under the condition that hBN with ${\gamma }_p=2{\mathrm{cm}}^{-1}$ is etched into a crystal with the same pattern as the metasurface. In this case, Rabi splitting is also observable, and strong coupling happens for the hBN crystal at any position, as indicated in Fig. 5(e). However, for the metasurface with hBN film or crystal at a specific position, the former always possesses a larger Rabi splitting and strong-coupling factor than the latter. This is because the hBN crystal has a smaller area than the film, so fewer phonons will get involved in the photon–phonon interaction, leading to the weaker strong coupling.

One of the most attractive features of topological photonic systems is their immunity to perturbations [8–10]. To investigate such topological protection characteristics of this metasurface, one square defect with side length $l^{\prime }$ is introduced in different positions of the supercell for every $3\times 3$ arrays, as presented in Fig. 6(a). Taking defect 1 as an example for the metasurface with no hBN, the transmission spectra under different defect lengths $l^{\prime }$ are plotted in Fig. 6(b). It is shown that as $l^{\prime }$ varies from 0 to 0.69 μm, the spectrum maintains a consistent line shape, and the transmission dip has a slight blueshift, introducing the frequency detuning $\mathrm{\Delta }$. In Fig. 6(c), we present the frequency detuning induced by defects in different positions. It is shown that the frequency shift will be more significant with a larger defect, but the maximum detuning is only $0.73{\mathrm{cm}}^{-1}$. The exception is defect 4, which has a negligible influence on the resonant frequency. These defects are further introduced in the hybrid metasurface. Here we only consider the largest defect because it will lead to the strongest perturbation in the strong coupling. With defects 1 to 3, the two emerged dips also exhibit a slight blueshift, as shown in Fig. 6(d); nevertheless, the Rabi splitting can still be identified, leading to the strong-coupling factor of 5.86, 6.45, and 6.15, respectively. Defect 4 causes little disturbance to strong coupling. Thus, the strong-coupling condition is still fulfilled in all cases. These results indicate that the corner state is topologically protected to a certain extent, and the induced tolerance further ensures the robustness of the strong coupling.

In the last section, we show that the resonant frequency of the corner state and the strong coupling can be dynamically tuned by stacking monolayer graphene film onto the metasurface. The surface conductivity model of graphene in the mid-infrared is described by the Drude-like model [73], $${\sigma }_{\mathrm{gra}}(\omega )=\frac{i{e}^2|E_F|}{\pi {\hslash }^2(\omega +i{\gamma }_e)}.$$(8)

Here, $e$ represents the electron charge, $E_F$ represents the Fermi level, $\hslash$ is the reduced Planck’s constant, and ${\gamma }_e$ is the electron damping rate set as $10{\mathrm{ps}}^{-1}$, according to previous reports [74]. The height of the metasurface $H$ is maintained at 4.245 μm.

The transmission spectra of the hybrid metasurface shown in the inset of Fig. 7(a) are plotted in Fig. 7(a). It is shown that the resonant frequency of the corner state remains at $1387.21{\mathrm{cm}}^{-1}$ when the Fermi level of graphene is at 0.1 eV; nonetheless, the $Q$ factor of the state is reduced from 3022 to 1950, compared with the case of no graphene, as indicated by the blue and red curves in Fig. 7(a). It is because graphene introduces additional losses that suppress the lifetime of the corner state. However, if the Fermi level is lifted to 1.0 eV, the resonant frequency will shift to $1388.48{\mathrm{cm}}^{-1}$, as indicated by the yellow curve in Fig. 7(a). Further blueshift of the resonant frequency happens when the Fermi level increases to 2.0 eV, as indicated by the purple curve in Fig. 7(a). We also notice that the $Q$ factor of the corner state will drop as the Fermi level increases because the higher Fermi level will lead to the larger real part of the graphene conductivity, as indicated in Eq. (8), which means more losses will be introduced. It should be mentioned that the Fermi level of the graphene can be flexibly tuned by applying voltage externally [41].

Further investigation shows that the resonant frequency of the corner state can be continuously tuned by changing the Fermi level of graphene, as presented in Fig. 7(b). With the range of the Fermi level from 0.2 to 2.0 eV, the resonant frequency can be continuously altered from 1387.56 to $1389.80{\mathrm{cm}}^{-1}$, with a significant frequency detuning of $2.24{\mathrm{cm}}^{-1}$. Such frequency variation can be attributed to the change of the effective permittivity of the graphene, which is caused by different Fermi levels. It will partly modify the permittivity of the surrounding environment and further influence the resonant frequency of the corner state, since the electric field of the corner state decays exponentially on the surface of the metasurface along the $z$ direction. According to simplified perturbation theory, the frequency detuning $\mathrm{\Delta }\omega$ of the corner state can be described as follows [75]: $$\frac{\mathrm{\Delta }\omega }{{\omega }_c}\propto \frac{\mathrm{\Delta }{\epsilon }_{\mathrm{gra}}{t}_{\mathrm{gra}}}{d_{\mathrm{field}}},$$(9)where $t_{\mathrm{gra}}$ is the equivalent thickness of the graphene, and $d_{\mathrm{field}}$ is the penetration depth of the electric field. ${\epsilon }_{\mathrm{gra}}$ is the in-plane effective permittivity of the graphene, which can be expressed as follows [76]: $${\epsilon }_{\mathrm{gra}}=2.5+\frac{i{\sigma }_{\mathrm{gra}}}{{\epsilon }_0\omega t_{\mathrm{gra}}}.$$(10)

Now we see that as the Fermi level rises, the real part of ${\epsilon }_{\mathrm{gra}}$ will also increase, resulting in positive $\mathrm{\Delta }{\epsilon }_{\mathrm{gra}}$ as well as the positive $\mathrm{\Delta }\omega$, which corresponds to the blueshift in the spectrum.

To gain more insights into this structure, we study the influence of graphene position along the $z$-axis on the frequency detuning. The relations between the resonant frequency of the corner state and Fermi level under different positions of graphene are presented in Fig. 7(c). It is shown that as the graphene comes from the top to the middle of the metasurface, the frequency detuning is gradually enhanced with the Fermi level from 0.2 to 2.0 eV. Especially, when graphene is in the middle, the frequency detuning is $6.36{\mathrm{cm}}^{-1}$, almost 3 times that when graphene is on the top. In addition, we also investigate the case when graphene crystal is on the top of the metasurface, and the consequent detuning is $2.05{\mathrm{cm}}^{-1}$, slightly lower than the case with graphene film on the top. All the results indicate that the larger overlap between graphene and the field of the corner state plays an important role in enhancing detuning performance.

Finally, the reconfigurable strong coupling is realized by assembling hBN (${\gamma }_p=2{\mathrm{cm}}^{-1}$) and graphene into the metasurface simultaneously, as illustrated in the inset of Fig. 7(d). Since the field of the corner state is symmetric to the $x-y$ plane, placing graphene at the bottom of the metasurface will exhibit the same tunability as placing it on the top. We tune the height $H$ to be 4.275 μm so that the resonant frequency of the corner state will match the phonon frequency when the Fermi level of the graphene is 1.0 eV. On this basis, the anticrossing behavior can be observed by altering the Fermi level from 0.2 to 2.0 eV, as displayed in Fig. 7(d), which arises from the frequency detuning of the corner state shown in Fig. 7(b). The strong-coupling factor here is 5.08, which meets the strong-coupling condition but is slightly lower than the case with no graphene. This attributes to the extra losses introduced by graphene, which increases the initial damping rate of the corner state. Such an hBN–metasurface–graphene heterostructure provides a feasible platform to realize the active control of strong coupling.

In summary, we propose a specially designed supercell metasurface composed of periodically arranged SOTIs based on a 2D SSH model, which supports nondegenerate topological multipolar corner states in the mid-infrared, and the two dipolar corner states, as well as their superposition state, can be directly excited by a far-field source due to the symmetry matching. By covering the metasurface with the hBN, topological phonon polaritons are realized through strong coupling between the corner state and phonon, indicated by Rabi splitting and anticrossing behavior. We further study the robustness of the corner state and strong coupling that originates from the topological protection. Finally, graphene is employed to realize dynamically reconfigurable corner states and strong coupling with on-demand electric control.

Extending isolated SOTI to be an array can promote novel large-area applications with topological protection. For example, with the assistance of the gain media, these collectively excited corner states can act as a large number of coherently coupled photonic emitters and form a laser array with large out-of-plane-emitting conformal aperture, which are of special importance for next-generation applications such as optical communication and light detection and ranging (lidar) [77,78]. Meanwhile, the thin analyte or nonlinear materials supported on the large-area surface of the array can interact sufficiently with light due to the in-plane field confinement, which can benefit the applications in label-free biosensing and optical nonlinear enhancement [79,80]. Moreover, recently emerging works have unveiled the connection between corner states and BICs [81], and nonlinear control of them has also been proposed [82]. Thus, the possibility to turn the two BICs in this metasurface system (monopole and quadrupole corner states) into leaky quasi-BICs is also worth exploring. Finally, the hybrid silicon-2D material metasurfaces have potential applications for topological polaritonic devices and reconfigurable topological devices in the mid-infrared region.