Since the introduction of the temporal coupled-mode theory (TCMT) into plasmonic optics by Fan^[1,2], vigorous developments of it have been witnessed. One of TCMT’s development paths is to simplify the incident light source (normally incident light) and detection (same polarization) but increase the number of modes and ports^[2-4], e.g. one-mode one-port TCMT^[5], one-mode two-port TCMT, two-mode one-port TCMT, two-mode two-port TCMT^[6,7]. TCMT has been used to explain a variety of novel phenomena, including Fano resonance^[8], plasmon-induced transparency (PIT)^[9], exceptional point (EP)^[10-14], and bound-state in the continuum^[15,16]. These phenomena have been widely used in coherent perfect absorbers^[17], unidirectional light propagation devices^[10,18], asymmetric mode switching^[19], sensor^[20-23], and improved Sagnac effect^[24]. Another development path is to keep the mode simple, while consider complicated source and detection (port) like oblique incident light sources^[25,26] and detection with different polarizations (polarization conversion)^[27]. In recent years, TCMT has been further generalized^[28-34]. For example, in 2021 Alu’s research group used TCMT to explain the thermal metasurface of the combined local and nonlocal light-matter interaction^[25]. In 2022, Fan’s research group used it to explain the optical force spectrum of resonant structure well^[32].

Particularly, the TCMT has been generalized to describe the cross-coupling terms of the near-field, which can be real^[6], virtual^[35,36], complex^[37], or chiral^[38]. The coupling can occur either between two bright modes^[31], or between the bright mode and dark mode^[6]. As revealed by TCMT, especially the two-mode two-port one, the total radiation and absorption losses are usually conserved irrespective to the change of structural asymmetry degree. In this paper, we studied the optical properties of the metasurfaces composed of a nano-bar and a nano-ring in a unit cell, where the total radiation loss of the system is found to be not conserved with varied asymmetry. There exist two degenerate points in the system. The mode corresponding to one degenerate point shows obvious localization, while the modes near the other degenerate point show obvious nonlocality. With the change of period, cross-coupling between the two modes is clearly shown, which is well explained by the generalized TCMT. This research provides a theoretical scheme for the new phenomenon of mode coupling and has certain guiding significance in the design of micro- and nano-photonics devices.

A unit cell of the metasurface is composed of a gold nano-ring and a gold nano-bar, as shown in Fig. 1. The bar with a thickness of $h_m=30\mathrm{nm}$ is inserted into the silica film with a thickness of $h=90\mathrm{nm}$, and the thickness of the ring is $h_m=30\mathrm{nm}$. Without losing generality, the period $P$ of the lattice is chosen to be 1000 nm. The length and width of the bar are $l=320\mathrm{nm}$ and $w=70\mathrm{nm}$, respectively. The outer and inner radii of the ring are $R=370\mathrm{nm}$ and $r=300\mathrm{nm}$, respectively. The asymmetry degree $d$ is defined as the lateral displacement of the bottom bar along the $y$-axis. In order to study the optical properties of metasurfaces, finite-difference time-domain (FDTD) simulations are performed. In the simulation, the refractive indices of the glass substrate and the intermediate layer are all set to $n=1.5$, and that for the background is 1. The dielectric constant of gold follows the Drude model $\epsilon =1-{\omega }_p^2/({\omega }^2+i\omega \gamma )$, where the frequencies of plasmonic resonances and damping are ${\omega }_p=1.374\times {10}^{16}\mathrm{rad}/\mathrm{s}$ and $\gamma =1.224\times {10}^{14}\mathrm{rad}/\mathrm{s}$.

It was shown previously^[39] that the PIT of this metasurface is caused by the coherent superposition of the electric dipole of the bar and the electric quadrupole of the ring, where the electric dipole acts as the bright mode and the electric quadrupole as the dark mode. The coupling between the two modes leads to the abnormal phenomenon of “merging-splitting-separation-proximity-merging”, i.e. the metasurfaces with varied $d$ present two frequency degeneracy points. One is EP, which refers to degeneracies of non-Hermitian systems with coalescent eigenstates and is commonly seen in systems with gain and loss such as parity-time-symmetric systems^[11]. The other one is the diabolic point (DP)^[40], which is the accidental degeneracy of energy levels in parameter space. What are the different optical properties near these two degenerate points? What causes the degeneracy point?

To answer the above questions, we use the FDTD approach to simulate spectra of metasurfaces with varied $d$. Figure 2 shows the simulated zeroth-order transmission ($T$), reflection ($R$), and absorption ($A$) spectra. The absorbance $A$ is calculated according to $$A=1-T-R.$$(1)

It can be seen from Fig. 2 that with the increased $d$, the transmission dip shows the “merging-splitting-separation-proximity-merging” and obvious PIT. In Fig. 2(c), A represents the absorption peak ($f_{\mathrm{A}}=180.4\mathrm{THz}$) of the metasurface with $d=60\mathrm{nm}$, whose $z$-component of the magnetic field in the unit cell on the $x=0\mathrm{nm}$ plane is shown in Fig. 2(d). The magnetic field is located around the rod, implying a local mode. B and C represent the two absorption peaks ($f_{\mathrm{B}}=189.1\mathrm{THz}$ and $f_{\mathrm{C}}=167.2\mathrm{THz}$) of the metasurface when $d=180\mathrm{nm}$, with the $z$-component distribution of the magnetic field shown in Figs. 2(e) and 2(f). The magnetic field of mode B spans the whole cell, and the exchange changes of field maximum and field minimum are shown along the $y$-direction. It is similar to the field distribution of surface plasmon polariton (SPP)^[41], therefore a global distribution. The magnetic field of mode C is also mainly localized near the bar, i.e., a local mode. D represents the absorption peak ($f_{\mathrm{D}}=173.7\mathrm{THz}$) (DP) of the metasurface with $d=360\mathrm{nm}$. The field distribution of the magnetic field, as plotted in Fig. 2(g), also shows a global field distribution similar to SPP. Why does the variation in $d$ change the properties of the mode? This is mainly because the distance between the bar and the ring becomes farther with $d=360\mathrm{nm}$. Once the electric dipole of the bar is excited and the dipolar field propagates and reaches the ring, the electric quadrupole of the ring is excited. Consequently, the field diffuses into the whole unit cell, showing the global characteristics. For $d=180\mathrm{nm}$, the long-wavelength C mode is the local mode and the short wavelength B mode is the global mode.

It can be seen from the above simulation that with changing $d$, the structure has two degenerate points. Are the optical properties of the structure near the two degenerate points the same? In order to answer this question, we first simulate the metasurfaces with asymmetry of $d=0$, 180, and 420 nm. The transmission, reflection, and absorption dispersion diagrams with the incident angle of plane wave varying from 0° to 20° are shown in Fig. 3. It can be seen from the figure that the mode with $d=0\mathrm{nm}$ is a local mode, and the resonance position does not shift significantly with the increase of incident angle. Two modes appear in the dispersion diagram with asymmetry $d=180\mathrm{nm}$, among which the long-wavelength one does not change obviously with the change of incident angle (local mode) and the short-wavelength one shows obvious redshift (global mode). At last, the two modes gradually merge into one single mode, which is particularly obvious in the absorption spectra. Moreover, with the increase of angle, the absorption increases significantly. In the dispersion diagram of $d=420\mathrm{nm}$, one mode is split into two modes. As the incident angle increases, the long-wavelength mode redshifts and the short-wavelength mode blueshifts, which show two obvious global modes. These are consistent with the conclusions given by the field distribution.

In order to further study the differences between the two frequency degeneracy points, we also simulated the periodic spectra of metasurfaces with $d=60\mathrm{nm}$ and $d=355\mathrm{nm}$ (other structural parameters are the same as above), as shown in Fig. 4. For $d=60\mathrm{nm}$, with the increase of the period, the quality factors of the two modes become larger and redshifts are observed. According to the Babinet complementarity principle^[42], SPP mode exists at a shorter wavelength, which is directly proportional to the change of period $P$; therefore, the short-wavelength mode is squeezed by the SPP at the shorter wavelength, resulting in a larger quality factor and redshift. The long-wavelength mode has a weak redshift under extrusion, and disappears due to the competition between two modes. For $d=355\mathrm{nm}$, when $P$ changes in steps of 20 nm from 920 nm to 1160 nm, the two modes first split, then merges into one mode, and then split again. To clearly see the evolution trend, we also plotted the transmission and reflection curves with $P$, as shown in Figs. 5(a) and 5(b). In the figure, we can see two modes appear splitting, merging, and then splitting with the increase of $P$. To explore whether this coupling is cross-coupling or anti cross-coupling, the field distributions of four modes with asymmetry $d=355\mathrm{nm}$ and period $P=1160\mathrm{nm}$ (or 920 nm) at wavelength ${\lambda }_{\mathrm{E}}=1.72\mathrm{\mu m}$ and ${\lambda }_{\mathrm{F}}=1.79\mathrm{\mu m}$ (or ${\lambda }_{\mathrm{G}}=1.63\mathrm{\mu m}$ and ${\lambda }_{\mathrm{H}}=1.73\mathrm{\mu m}$) are simulated. Figures 5(c) and 5(d) show the field distribution of the asymmetric electric quadrupole in mode E, and the field in the lower half of the ring is relatively strong. Mode F is also the field of the electric quadrupole, but it is mainly concentrated in the upper half ring. Figures 5(e) and 5(f) show the field distribution of modes G and H, also corresponding to the asymmetric electric quadrupole, but with the field concentration opposite to those of the modes E and F. That is, with the change of period, mode E evolves into mode H and mode F evolves into mode G, and the two modes appear cross-coupling. The reason is that DP is the degeneracy point between modes due to the complex and subtle relationship between the period and the size of metasurfaces.

The results show that the mode of the structures near EP related to the radiation loss of the bright mode is almost nondispersive. With the increase of $P$, the mode has a weak redshift, and the quality factor increases due to the extrusion caused by the redshifts of other modes. The DP due to the adjustment of structural parameters has obvious redshifts with the change of angle. With the increase of $P$, the phenomenon of splitting-merging-splitting occurs, and cross-coupling appears. The mode E becomes the mode H and the mode F becomes the mode G. According to Ref. [6], these phenomena can be explained by TCMT. The dipole of the bar acts as the bright mode (p) and the quadrupole of the ring acts as the dark mode (q). With the symmetry breaking of the structures, the coherent superposition of the two modes leads to the asymmetric line shape of the spectra. To further study the physical mechanism behind it, we choose the two-mode two-port TCMT. The Hamiltonian equation without coupling between two modes can be written as $${\widehat{\mathrm{H}}}_0|{\mathrm{\Psi }}_{\mathrm{p}}\rangle =2\pi f_{\mathrm{p}}|{\mathrm{\Psi }}_{\mathrm{p}}\rangle ,$$(2)$${\widehat{\mathrm{H}}}_0|{\mathrm{\Psi }}_{\mathrm{q}}\rangle =2\pi f_{\mathrm{q}}|{\mathrm{\Psi }}_{\mathrm{q}}\rangle ,$$(3)where $|{\mathrm{\Psi }}_{\mathrm{p}}\rangle$ and $|{\mathrm{\Psi }}_{\mathrm{q}}\rangle$ are eigenstates of the dipole and quadrupole modes for the symmetric nanostructure without coupling. $f_{\mathrm{p}}$ and $f_{\mathrm{q}}$ are their corresponding frequency. Consider our plasmonic resonator as a two-mode two-port system. And we define $S_j^{+}$ and $S_j^{-}$ ($j=1,2$) are the complex amplitude of the incoming and outgoing waves at the $j$th port. The amplitudes of the dipole and quadrupole modes in our nanostructure are represented by $a_{\mathrm{p}}$ and $a_{\mathrm{q}}$, respectively. The dynamic equations for the amplitude of the resonance mode can be written as $$\begin{gathered} \frac{1}{2\pi }\frac{\mathrm{d}}{\mathrm{d}t}\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right)=-i[\left(\begin{array}{cc}{f}_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}& {\kappa }_{\mathrm{p}\mathrm{q}}\\ {\kappa }_{\mathrm{q}\mathrm{p}}& f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}\end{array}\right)-i\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}^{\prime }& 0\\ 0& {\gamma }_{\mathrm{q}}^{\prime }\end{array}\right)]\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right) \\ -\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}& 0\\ 0& 0\end{array}\right)\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right)+\left(\begin{array}{cc}{k}_{\mathrm{p}1}& k_{\mathrm{p}2}\\ 0& 0\end{array}\right)\left(\begin{array}{c}{S}_1^{+}\\ S_2^{+}\end{array}\right), \end{gathered}$$(4)$$\left(\begin{array}{c}{S}_1^{-}\\ S_2^{-}\end{array}\right)=\left(\begin{array}{cc}{S}_{11}^0& S_{12}^0\\ S_{21}^0& S_{22}^0\end{array}\right)\left(\begin{array}{c}{S}_1^{+}\\ S_2^{+}\end{array}\right)+\left(\begin{array}{cc}{d}_{1p}& 0\\ d_{2p}& 0\end{array}\right)\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right),$$(5)where ${\kappa }_{\mathrm{p}\mathrm{p}}$ and ${\kappa }_{\mathrm{q}\mathrm{q}}$ denote the on-site corrections by the perturbation and ${\kappa }_{\mathrm{p}\mathrm{q}}={\kappa }_{\mathrm{q}\mathrm{p}}$ describes the near-field coupling between two original modes, again contributed by the perturbation. ${\gamma }_{\beta }$ (or ${\gamma }_{\beta }^{\prime }$) is the radiation (or absorption) loss rate for the $\beta$th mode. $k_{\beta j}$ or $d_{j\beta }(j=1,2)$ is the coupling constant between the $\beta$th mode and the incoming or outgoing wave at the $j$th port. $S_{jl}^0(j,l=1,2)$ is the $S$-parameter matrix element of the background system. The coupling constant $k_{\beta j}$ or $d_{j\beta }$ and the radiation loss rate ${\gamma }_{\beta }$ are dependent on each other by energy conservation and time-reversal symmetry consideration. According to principles of optics^[43], under normal incidence or grazing incidence, the half-wave loss will occur when light is reflected from light sparse to light dense medium. Let medium 2 be glass and medium 1 be air. The $S$-parameter matrix elements of the background system are $$S_{12}^0=S_{21}^0=\frac{2\sqrt{n}}{n+1},S_{22}^0=-S_{11}^0=\frac{n-1}{n+1}.$$(6)

We diagonalize the near-field term to obtain a unitary matrix and then carry out a unitary transformation for Eqs. (4) and (5). At first, the eigenfrequencies are obtained: $${\tilde{f}}_{\mp }=\overline{f}\mp \mathrm{\Delta }\tilde{f}.$$(7)

Here, we have defined $\overline{f}=(f_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}+f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}})/2$, $\mathrm{\Delta }\tilde{f}=\sqrt{{(\mathrm{\Delta }f)}^2+{\kappa }_{\mathrm{p}\mathrm{q}}^2}$ and $\mathrm{\Delta }f=(f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}-f_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}})/2$ as the frequency average and difference between two modes; the subscripts “$+$” and “$-$” correspond to the “dressed” dark and bright modes, respectively. We can obtain the eigenvector of the “dressed” dark mode and the “dressed” bright mode, and the transformation matrix as follows: $$M=\left(\begin{array}{cc}\frac{-\mathrm{\Delta }f-\mathrm{\Delta }\tilde{f}}{\sqrt{{(-\mathrm{\Delta }f-\mathrm{\Delta }\tilde{f})}^2+{\kappa }_{\mathrm{p}\mathrm{q}}^2}}& \frac{-\mathrm{\Delta }f+\mathrm{\Delta }\tilde{f}}{\sqrt{{(-\mathrm{\Delta }f+\mathrm{\Delta }\tilde{f})}^2+{\kappa }_{\mathrm{p}\mathrm{q}}^2}}\\ \frac{{\kappa }_{\mathrm{p}\mathrm{q}}}{\sqrt{{(-\mathrm{\Delta }f-\mathrm{\Delta }\tilde{f})}^2+{\kappa }_{\mathrm{p}\mathrm{q}}^2}}& \frac{{\kappa }_{\mathrm{p}\mathrm{q}}}{\sqrt{{(-\mathrm{\Delta }f+\mathrm{\Delta }\tilde{f})}^2+{\kappa }_{\mathrm{p}\mathrm{q}}^2}}\end{array}\right).$$(8)

Through unitary transformations on both sides of Eqs. (4) and (5), we obtain the dynamical equation of coupled representation as follows: $$\begin{gathered} \frac{1}{2\pi }\frac{\mathrm{d}}{\mathrm{d}t}\left(\begin{array}{c}{\tilde{a}}_{+}\\ {\tilde{a}}_{-}\end{array}\right)=-i\left(\begin{array}{cc}{\tilde{f}}_{+}& 0\\ 0& {\tilde{f}}_{-}\end{array}\right)\left(\begin{array}{c}{\tilde{a}}_{+}\\ {\tilde{a}}_{-}\end{array}\right) \\ +\left(\begin{array}{cc}{\tilde{\gamma }}_{+}+{{\tilde{\gamma }}^{\prime }}_{+}& \tilde{X}\\ \tilde{X}& {\tilde{\gamma }}_{-}+{{\tilde{\gamma }}^{\prime }}_{-}\end{array}\right)\left(\begin{array}{c}{\tilde{a}}_{+}\\ {\tilde{a}}_{-}\end{array}\right)+\left(\begin{array}{cc}{\tilde{k}}_{+1}& {\tilde{k}}_{+2}\\ {\tilde{k}}_{-1}& {\tilde{k}}_{-2}\end{array}\right)\left(\begin{array}{c}{S}_{10}^{+}\\ S_{20}^{+}\end{array}\right). \end{gathered}$$(9)

The output wave equation of coupled representation is$$\left(\begin{array}{c}{S}_{10}^{-}\\ S_{20}^{-}\end{array}\right)=\left(\begin{array}{cc}{S}_{11}^0& S_{12}^0\\ S_{21}^0& S_{22}^0\end{array}\right)\left(\begin{array}{c}{S}_{10}^{+}\\ S_{20}^{+}\end{array}\right)+\left(\begin{array}{cc}{\tilde{d}}_{1+}& {\tilde{d}}_{1-}\\ {\tilde{d}}_{2+}& {\tilde{d}}_{2+}\end{array}\right)\left(\begin{array}{c}{\tilde{a}}_{+}\\ {\tilde{a}}_{-}\end{array}\right).$$(10)

The expressions of physical quantities in coupled representation are obtained: $$\tilde{X}=-{\kappa }_{\mathrm{p}\mathrm{q}}{\gamma }_{\mathrm{p}}/(2\mathrm{\Delta }\tilde{f}),$$(11)$${\tilde{\gamma }}_{-}={\gamma }_{\mathrm{p}}(\mathrm{\Delta }\tilde{f}+\mathrm{\Delta }f)/2\mathrm{\Delta }\tilde{f},$$(12)$${\tilde{\gamma }}_{+}={\gamma }_{\mathrm{p}}(\mathrm{\Delta }\tilde{f}-\mathrm{\Delta }f)/2\mathrm{\Delta }\tilde{f},$$(13)$${\tilde{\gamma }}_{-}^{\prime }=\frac{{\gamma }_{\mathrm{p}}^{\prime }+{\gamma }_{\mathrm{q}}^{\prime }}{2}-\frac{({\gamma }_{\mathrm{q}}^{\prime }-{\gamma }_{\mathrm{p}}^{\prime })\mathrm{\Delta }f}{2\mathrm{\Delta }\tilde{f}},$$(14)$${{\tilde{\gamma }}^{\prime }}_{+}=\frac{{\gamma }_{\mathrm{p}}^{\prime }+{\gamma }_{\mathrm{q}}^{\prime }}{2}+\frac{({\gamma }_{\mathrm{q}}^{\prime }-{\gamma }_{\mathrm{p}}^{\prime })\mathrm{\Delta }f}{2\mathrm{\Delta }\tilde{f}}.$$(15)

According to Eqs. (12) to (15), we obtain two conserved quantities: $${\tilde{\gamma }}_{+}+{\tilde{\gamma }}_{-}={\gamma }_{\mathrm{p}},$$(16)$${{\tilde{\gamma }}^{\prime }}_{+}+{{\tilde{\gamma }}^{\prime }}_{-}={\gamma }_{\mathrm{p}}^{\prime }+{\gamma }_{\mathrm{q}}^{\prime }.$$(17)

Equations (16) and (17) show the total radiation loss and the total absorption loss are both conserved quantities during the whole evolution of the system. At the same time, Eqs. (7), (12), and (13) can inversely solve the near-field coupling coefficient under the uncoupled basis vector. $${\kappa }_{\mathrm{p}\mathrm{q}}=\frac{{\tilde{f}}_{+}-{\tilde{f}}_{-}}{2{\gamma }_{\mathrm{p}}}\sqrt{({\gamma }_{\mathrm{p}}+{\tilde{\gamma }}_{+}-{\tilde{\gamma }}_{-})({\gamma }_{\mathrm{p}}-{\tilde{\gamma }}_{+}+{\tilde{\gamma }}_{-})},$$(18)$${\kappa }_{\mathrm{q}\mathrm{q}}=[{\tilde{f}}_{+}+{\tilde{f}}_{-}-f_{\mathrm{p}}-f_{\mathrm{q}}+\frac{{\tilde{f}}_{+}-{\tilde{f}}_{-}}{{\gamma }_{\mathrm{p}}}({\tilde{\gamma }}_{-}-{\tilde{\gamma }}_{+})]/2,$$(19)$${\kappa }_{\mathrm{p}\mathrm{p}}=[{\tilde{f}}_{+}+{\tilde{f}}_{-}-f_{\mathrm{p}}-f_{\mathrm{q}}-\frac{{\tilde{f}}_{+}-{\tilde{f}}_{-}}{{\gamma }_{\mathrm{p}}}({\tilde{\gamma }}_{-}-{\tilde{\gamma }}_{+})]/2.$$(20)

The complex amplitude of transmission under the coupling basis vector can be obtained by $$t=\frac{2\sqrt{n}}{n+1}-\frac{2\sqrt{n}}{n+1}(\frac{{\tilde{\gamma }}_{-}}{{\tilde{W}}_{-}}+\frac{{\tilde{\gamma }}_{+}}{{\tilde{W}}_{+}}-2\frac{{\tilde{\gamma }}_{-}{\tilde{\gamma }}_{+}}{{\tilde{W}}_{-}{\tilde{W}}_{+}})/(1-\frac{{\tilde{\gamma }}_{-}{\tilde{\gamma }}_{+}}{{\tilde{W}}_{-}{\tilde{W}}_{+}}),$$(21)where ${\tilde{W}}_{\mp }=[-i(f-{\tilde{f}}_{\mp })+{\tilde{\gamma }}_{\mp }+{\tilde{\gamma }}_{\mp }^{\prime }]$. Based on Eq. (21), the simulation results are fitted. The red circle represents the simulation results and the blue solid curve for the fitting ones, as shown in Fig. 6. The results show that the theory is in good agreement with the simulation. Figure 7 is the variation curve of fitting parameters with $d$. Figure 7(a) is the variation of the resonance frequency of two dressed modes with $d$, showing the phenomenon of “merging-splitting-separation-proximity-merging”. Figure 7(b) shows that the radiation loss of the two dressed modes changes with $d$. We can see the radiation loss of the dressed bright mode first decreases rapidly from the maximum and then decreases slowly, and that of the dressed dark mode increases from 0 and then decreases slowly, showing a phenomenon different from the common increase. Figure 7(c) is the curve of the absorption loss of the two dressed modes with $d$. It can be seen from the figure that the absorption of the two modes is almost unchanged. Figure 7(d) shows the sum frequency ${\tilde{f}}_t$ of the two dressed modes varies with $d$. It can be seen from the figure that ${\tilde{f}}_t$ decreases slowly before $d<360\mathrm{nm}$, decreases suddenly at 360 nm, and then decreases slowly. Figure 7(e) shows the total radiation losses ${\tilde{\gamma }}_t$ with $d$. It is obvious from the figure that when $d<200\mathrm{nm}$, ${\tilde{\gamma }}_t$ is almost conserved, about 10.3 THz. When $d>200\mathrm{nm}$, ${\tilde{\gamma }}_t$ decreases linearly, and finally to the minimum of 3.7 THz, as shown in the shadow in Fig. 7(e). That is, when $d>200\mathrm{nm}$, the total radiation loss is not conserved, which is inconsistent with the conservation of the total radiation loss when a near-field coupling coefficient is a real number. Figure 7(f) shows the total absorption loss ${\tilde{\gamma }}_t^{\prime }$ with $d$. It can be seen from the figure that ${\tilde{\gamma }}_t^{\prime }$ remains almost constant, which is consistent with the conservation of the total absorption loss when a near-field coupling coefficient is a real number.

Why is the total radiation loss not conserved when $d>200\mathrm{nm}$? How to explain it? When $d>200\mathrm{nm}$, the bar and the ring are far away. The dipole of the bar is excited and propagates a certain distance to the ring. The quadrupole of the ring is excited by near-field coupling. Because the ring is far from the bar, there is a propagation phase in the near-field cross-coupling coefficient. The on-site term of near-field coupling indicates that the field of the bar propagates to the ring, and then the quadrupole of the ring is excited. The quadrupole field of the ring has a certain perturbation to the dipole field of the bar. Therefore, the propagation phase corresponding to this item is larger and it is reasonable to introduce the propagation phase into the near-field coupling coefficient. The master equation of dynamics under the uncoupling basis vector is $$\begin{gathered} \frac{1}{2\pi }\frac{\mathrm{d}}{\mathrm{d}t}\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right)=-i[\left(\begin{array}{cc}{f}_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}{e}^{-i\alpha }& {\kappa }_{\mathrm{p}\mathrm{q}}{e}^{-i\beta }\\ {\kappa }_{\mathrm{p}\mathrm{q}}{e}^{-i\beta }& f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}{e}^{-i\alpha }\end{array}\right)-i\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}^{\prime }& 0\\ 0& {\gamma }_{\mathrm{q}}^{\prime }\end{array}\right)]\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right) \\ -\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}& 0\\ 0& 0\end{array}\right)\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right)+\left(\begin{array}{cc}{k}_{\mathrm{p}1}& k_{\mathrm{p}2}\\ k_{\mathrm{q}1}& k_{\mathrm{q}2}\end{array}\right)\left(\begin{array}{c}{S}_{10}^{+}\\ S_{20}^{+}\end{array}\right) \\ =-i[\left(\begin{array}{cc}{f}_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\cos \alpha & {\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta \\ {\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta & f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}\cos \alpha \end{array}\right)-i\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}^{\prime }& 0\\ 0& {\gamma }_{\mathrm{q}}^{\prime }\end{array}\right)]\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right) \\ -\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha & {\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta \\ {\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta & {\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha \end{array}\right)\left(\begin{array}{c}{a}_{\mathrm{p}}\\ a_{\mathrm{q}}\end{array}\right)+\left(\begin{array}{cc}{k}_{\mathrm{p}1}& k_{\mathrm{p}2}\\ k_{\mathrm{q}1}& k_{\mathrm{q}2}\end{array}\right)\left(\begin{array}{c}{S}_{10}^{+}\\ S_{20}^{+}\end{array}\right), \end{gathered}$$(22)where $\left(\begin{array}{cc}{f}_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\cos \alpha & {\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta \\ {\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta & f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}\cos \alpha \end{array}\right)$ is the near-field coupling matrix, $\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha & {\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta \\ {\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta & {\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha \end{array}\right)$ is the radiation loss matrix, and $\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}^{\prime }& 0\\ 0& {\gamma }_{\mathrm{q}}^{\prime }\end{array}\right)$ is the absorption loss matrix. For simplicity, we set ${\kappa }_{\mathrm{p}\mathrm{p}}\cos \alpha ={\kappa }_{\mathrm{p}\mathrm{p}}^{\prime }$, ${\kappa }_{\mathrm{q}\mathrm{q}}\cos \alpha ={\kappa }_{\mathrm{q}\mathrm{q}}^{\prime }$, and ${\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta ={\kappa }_{\mathrm{p}\mathrm{q}}^{\prime }$. Because the coupling coefficient between the incident port and two modes is not zero, both modes are bright modes. Similarly, the coupling coefficient between two modes and the exit port is also not zero. The eigenvalues can be obtained by diagonalizing the near-field coupling term $${{\tilde{f}}^{\prime }}_{\mp }={\overline{f}}^{\prime }\mp \mathrm{\Delta }{\tilde{f}}^{\prime }.$$(23)

Here, we have defined ${\overline{f}}^{\prime }=(f_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}^{\prime }+f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}^{\prime })/2$, $\mathrm{\Delta }{\tilde{f}}^{\prime }=\sqrt{{(\mathrm{\Delta }{f}^{\prime })}^2+{{\kappa }^{\prime }}_{\mathrm{p}\mathrm{q}}^2}$ and $\mathrm{\Delta }{f}^{\prime }=(f_{\mathrm{q}}+{\kappa }_{\mathrm{q}\mathrm{q}}^{\prime }-f_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}}^{\prime })/2$ as the frequency average and difference between two modes, the subscript “$+$” corresponds to the dressed dark mode, and the subscript “$-$” corresponds to the dressed bright mode. We substitute Eq. (23) into eigenvalue equations, and we can obtain the eigenvectors of the dressed bright mode and the dressed dark mode, and then obtain the transformation matrix as follows: $${\tilde{M}}^{\prime }=\left(\begin{array}{cc}-\frac{\mathrm{\Delta }{f}^{\prime }+\mathrm{\Delta }{\tilde{f}}^{\prime }}{2{\kappa }_{\mathrm{p}\mathrm{q}}^{\prime }}& \frac{-\mathrm{\Delta }{f}^{\prime }+\mathrm{\Delta }{\tilde{f}}^{\prime }}{2{\kappa }_{\mathrm{p}\mathrm{q}}^{\prime }}\\ 1& 1\end{array}\right).$$(24)

Through unitary transformations, the expressions of physical quantities in coupled representation are obtained: $$\left(\begin{array}{cc}{\tilde{\gamma }}_{-}& \tilde{X}\\ \tilde{X}& {\tilde{\gamma }}_{+}\end{array}\right)\equiv {M^{\prime }}^{-1}\left(\begin{array}{cc}{\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha & {\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta \\ {\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta & {\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha \end{array}\right)M^{\prime }=\left(\begin{array}{cc}{B}_{11}& B_{12}\\ B_{21}& B_{22}\end{array}\right),$$(25)where $$B_{11}=\frac{\mathrm{\Delta }{f}^{\prime }({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha -{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )+\mathrm{\Delta }{\tilde{f}}^{\prime }({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha -2{\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )}{2\mathrm{\Delta }{\tilde{f}}^{\prime }},$$(26)$$B_{12}=\frac{{\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta (-{\gamma }_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )}{2\mathrm{\Delta }{\tilde{f}}^{\prime }},$$(27)$$B_{21}=\frac{{\kappa }_{\mathrm{p}\mathrm{q}}\cos \beta (-{\gamma }_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )}{\mathrm{\Delta }{\tilde{f}}^{\prime }}-\frac{2\mathrm{\Delta }{f}^{\prime }\sin \beta }{\cos \beta },$$(28)$$B_{22}=\frac{\mathrm{\Delta }{f}^{\prime }(-{\gamma }_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )+\mathrm{\Delta }{\tilde{f}}^{\prime }({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +2{\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )}{2\mathrm{\Delta }{\tilde{f}}^{\prime }}.$$(29)

In order to meet the requirements of energy conservation and time inversion symmetry in the coupled TCMT, here $B_{12}$ must be equal to $B_{21}$, then $$\frac{2\mathrm{\Delta }{f}^{\prime }\sin \beta }{\cos \beta }=0.$$(30)

Discuss: 1.When $\sin \beta =0$, we obtain expressions of physical quantities in coupled representation: $$\tilde{X}=-{\kappa }_{\mathrm{p}\mathrm{q}}({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha -{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )/(2\mathrm{\Delta }{\tilde{f}}^{\prime }),$$(31)$${\tilde{\gamma }}_{-}=\frac{\mathrm{\Delta }{f}^{\prime }({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha -{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )+\mathrm{\Delta }{\tilde{f}}^{\prime }({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )}{2\mathrm{\Delta }{\tilde{f}}^{\prime }},$$(32)$${\tilde{\gamma }}_{+}=\frac{\mathrm{\Delta }{f}^{\prime }(-{\gamma }_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )+\mathrm{\Delta }{\tilde{f}}^{\prime }({\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha )}{2\mathrm{\Delta }{\tilde{f}}^{\prime }}.$$(33)When $\alpha =0$, Eqs. (32) and (33) degenerate into Eqs. (12) and (13).2.When $\mathrm{\Delta }{f}^{\prime }=0$, then $$\tilde{X}=\frac{-{\gamma }_{\mathrm{p}}-{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha }{2},$$(34)$${\tilde{\gamma }}_{-}=\frac{{\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha -2{\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha }{2},$$(35)$${\tilde{\gamma }}_{+}=\frac{{\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +2{\kappa }_{\mathrm{p}\mathrm{q}}\sin \beta +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha }{2}.$$(36)

Generally speaking, it is difficult to ensure $\mathrm{\Delta }{f}^{\prime }=0$ during the whole evolution process, so here $\sin \beta =0$ is selected. The absorption losses under the coupling basis vector are $${{\tilde{\gamma }}^{\prime }}_{+}=\frac{({\gamma }_{\mathrm{p}}^{\prime }-{\gamma }_{\mathrm{q}}^{\prime })\mathrm{\Delta }f}{2\mathrm{\Delta }\tilde{f}}+\frac{{\gamma }_{\mathrm{p}}^{\prime }+{\gamma }_{\mathrm{q}}^{\prime }}{2},$$(37)$${{\tilde{\gamma }}^{\prime }}_{-}=\frac{{\gamma }_{\mathrm{p}}^{\prime }+{\gamma }_{\mathrm{q}}^{\prime }}{2}-\frac{({\gamma }_{\mathrm{p}}^{\prime }-{\gamma }_{\mathrm{q}}^{\prime })\mathrm{\Delta }f}{2\mathrm{\Delta }\tilde{f}}.$$(38)

Based on Eqs. (32), (33) and (37), (38), we obtain $${\tilde{\gamma }}_{+}+{\tilde{\gamma }}_{-}={\gamma }_{\mathrm{p}}+{\kappa }_{\mathrm{p}\mathrm{p}}\sin \alpha +{\kappa }_{\mathrm{q}\mathrm{q}}\sin \alpha ,$$(39)$${{\tilde{\gamma }}^{\prime }}_{+}+{{\tilde{\gamma }}^{\prime }}_{-}={\gamma }_{\mathrm{p}}^{\prime }+{\gamma }_{\mathrm{q}}^{\prime }.$$(40)

Equation (40) shows the total absorption loss is a conserved quantity during the whole evolution of the system. However, Eq. (39) shows the total radiation loss is not a conserved quantity. Combining Eqs. (23) and (32), (33), we obtain $$\tan \alpha =\frac{{\tilde{\gamma }}_{+}+{\tilde{\gamma }}_{-}-{\gamma }_{\mathrm{p}}}{{\tilde{f}}_{-}+{\tilde{f}}_{+}-f_{\mathrm{p}}-f_{\mathrm{q}}}.$$(41)

When $f_{\mathrm{q}}=f_{\mathrm{p}}$, further derivation can be obtained: $${\kappa }_{\mathrm{p}\mathrm{p}}=\frac{1}{2}[\frac{{\tilde{\gamma }}_{+}+{\tilde{\gamma }}_{-}-{\gamma }_{\mathrm{p}}}{\sin \alpha }-\frac{{\gamma }_{\mathrm{p}}\cos \alpha \mp \sqrt{{({\gamma }_{\mathrm{p}}\cos \alpha )}^2-4\sin \alpha \cos \alpha ({\tilde{\gamma }}_{-}-{\tilde{\gamma }}_{+})({\tilde{f}}_{+}-{\tilde{f}}_{-})}}{2\sin \alpha \cos \alpha }],$$(42)$${\kappa }_{\mathrm{q}\mathrm{q}}=\frac{1}{2}[\frac{{\tilde{\gamma }}_{+}+{\tilde{\gamma }}_{-}-{\gamma }_{\mathrm{p}}}{\sin \alpha }+\frac{{\gamma }_{\mathrm{p}}\cos \alpha \mp \sqrt{{({\gamma }_{\mathrm{p}}\cos \alpha )}^2-4\sin \alpha \cos \alpha ({\tilde{\gamma }}_{-}-{\tilde{\gamma }}_{+})({\tilde{f}}_{+}-{\tilde{f}}_{-})}}{2\sin \alpha \cos \alpha }],$$(43)$${\kappa }_{\mathrm{p}\mathrm{q}}=\frac{1}{2}{\left\{{({\tilde{f}}_{-}-{\tilde{f}}_{+})}^2-{\left[\frac{{\gamma }_{\mathrm{p}}\cos \alpha \mp \sqrt{{({\gamma }_{\mathrm{p}}\cos \alpha )}^2-4\sin \alpha \cos \alpha ({\tilde{\gamma }}_{-}-{\tilde{\gamma }}_{+})({\tilde{f}}_{-}-{\tilde{f}}_{+})}}{2\sin \alpha }\right]}^2\right\}}^{0.5}.$$(44)

When $\alpha =0$, Eqs. (42)–(44) degenerate to Eqs. (18)–(20). Therefore, our theory is self-consistent.

In order to verify the correctness of our hypothesis, based on Eqs. (18)–(20) [Eqs. (42)–(44)], we plot the near-field coupling coefficient of metasurfaces with $d\le 180\mathrm{nm}$ (or $d\ge 180\mathrm{nm}$), as shown in Fig. 8, where $f_{\mathrm{p}}=f_{\mathrm{q}}=179.2\mathrm{THz}$, ${\gamma }_{\mathrm{p}}=10.3\mathrm{THz}$, ${\gamma }_{\mathrm{p}}^{\prime }=3.6\mathrm{THz}$, and ${\gamma }_{\mathrm{q}}^{\prime }=3.5\mathrm{THz}$. Figures 8(a)–8(c) show the near-field coupling coefficients of metasurfaces with $d\le 180\mathrm{nm}$. As can be seen from the figure, ${\kappa }_{\mathrm{p}\mathrm{q}}$ increases linearly from 0 to 10.67 THz, ${\kappa }_{\mathrm{p}\mathrm{p}}$ decreases almost linearly from 0 to $-4.4\mathrm{THz}$, and ${\kappa }_{\mathrm{q}\mathrm{q}}$ slowly increases from 0 to 1.7 THz. Figures 8(d)–8(f) show the near-field coupling coefficient ${\kappa }_{\mathrm{p}\mathrm{q}}$, the real part of ${\kappa }_{\mathrm{p}\mathrm{p}}$, and the real part of ${\kappa }_{\mathrm{q}\mathrm{q}}$ for $d\ge 180\mathrm{nm}$, respectively. It can be seen from the figure that the near-field coupling coefficient ${\kappa }_{\mathrm{p}\mathrm{q}}$ with $d\ge 180\mathrm{nm}$ decreases almost linearly from 10.66 THz to 0 with $d=360\mathrm{nm}$. When $d\ge 360\mathrm{nm}$, ${\kappa }_{\mathrm{p}\mathrm{q}}$ remains 0. The real part of ${\kappa }_{\mathrm{p}\mathrm{p}}$ decreases linearly from $-4.4\mathrm{THz}$ to $-8.8\mathrm{THz}$, and here has small fluctuations around $-7\mathrm{THz}$. The real part of ${\kappa }_{\mathrm{q}\mathrm{q}}$ increases slowly from 1.8 THz to 3.5 THz, and then decreases to 0 with $d=340\mathrm{nm}$. When $d=360\mathrm{nm}$, ${\kappa }_{\mathrm{q}\mathrm{q}}$ changes suddenly to $-6.9\mathrm{THz}$, and then slowly decreases to $-8.8\mathrm{THz}$ with the increase of $d$. At $d=180\mathrm{nm}$, the difference between the two methods is very small, which also shows that the generalized theory is self-consistent.

Figure 9(a) shows the propagation phase $\alpha$ variation with $d$. As can be seen from the figure, when $d<360\mathrm{nm}$, $\alpha$ increases linearly from 0 to 0.77. When $d=360\mathrm{nm}$, $\alpha$ changes suddenly from 0.77 to 0.46, then slowly to 0.415. Figure 9(b) shows the imaginary part of $c\kappa ={\kappa }_{\mathrm{p}\mathrm{p}}-{\kappa }_{\mathrm{q}\mathrm{q}}$ with $d$. It can be seen from the figure that $\mathrm{Im}c\kappa$ when $240\mathrm{nm}<d<360\mathrm{nm}$ is not 0. Because the existence of the propagation phase makes the near-field coupling term more complex. The existence of the propagation phase not only changes the resonance frequency of dressed mode, but also changes the total radiation loss of dressed mode. Theoretical studies show that DP is related to the propagation phase, which is related to the period and the size of nanostructures. Therefore, the generalized TCMT well explains the simulation phenomenon.

In recent years, TCMT has made a great breakthrough, successfully explained various optical phenomena, and has a wide range of application prospects. In this paper, the optical properties of metasurfaces with the unit cell composed of nano-bar and nano-ring are numerically simulated. The TCMT under the coupled basis vector well explains the numerical simulation results. However, when the asymmetry degree is large, the total radiation loss is not conserved, which is contrary to the traditional TCMT requirement. The generalized TCMT by introducing the propagation phase into the near-field coupling term, solves this contradiction well. The results show that, unlike the local mode near EP, the global mode near DP is related to the propagation phase and shows the phenomenon of cross-coupling with the change of period. This study provides a new theoretical scheme for understanding the interaction between light and matter and has certain guiding significance for popularizing it in various related fields.