Intercoupling of Cascaded Metasurfaces for Broadband Spectral Scalability

Electromagnetic metasurfaces have been intensively used as ultra-compact and easy-to-integrate platforms for versatile wave manipulations from optical to terahertz (THz) and millimeter wave (MMW) ranges. In this paper, the less investigated effects of the interlayer coupling of multiple metasurfaces cascaded in parallel are intensively exploited and leveraged for scalable broadband spectral regulations. The hybridized resonant modes of cascaded metasurfaces with interlayer couplings are well interpreted and simply modeled by the transmission line lumped equivalent circuits, which are used in return to guide the design of the tunable spectral response. In particular, the interlayer gaps and other parameters of double or triple metasurfaces are deliberately leveraged to tune the inter-couplings for as-required spectral properties, i.e., the bandwidth scaling and central frequency shift. As a proof of concept, the scalable broadband transmissive spectra are demonstrated in the millimeter wave (MMW) range by cascading multilayers of metasurfaces sandwiched together in parallel with low-loss dielectrics (Rogers 3003). Finally, both the numerical and experimental results confirm the effectiveness of our cascaded model of multiple metasurfaces for broadband spectral tuning from a narrow band centered at 50 GHz to a broadened range of 40~55 GHz with ideal side steepness, respectively.


Analysis Model of Single-Layer Metasurface
To probe the inter-coupling effect between multiple metasurfaces cascaded in parallel, we start with the accurate modeling of a single metasurface. First, one type of metasurface is constructed by using the Jerusalem cross inside a square ring as the unit meta-atom cell ( Figure 1a). Typically, metasurfaces with such simple capacitive cells of metal elements act as frequency selective devices with transmissive band-stop spectra, as shown in Figure 1b, calculated by the finite element method (FEM). meta-atom cell ( Figure 1a). Typically, metasurfaces with such simple capacitive cells of metal elements act as frequency selective devices with transmissive band-stop spectra, as shown in Figure 1b, calculated by the finite element method (FEM). Constructed with symmetrically and periodically arranged elements, such a metasurface shows an ideal geometrical isotropy and polarization independence. So, two basic lumped equivalent circuits with similar band-stop transmission spectra can be used for simple modeling: the parallel and series configurations of LC circuits, as shown in Figure 2, according to the transmission line (TL) theory [51]. Essentially, these two circuits differ in the imaginary parts of their characteristic impedance. Namely, the dominant reactances of the parallel and series LC circuit are − /( − ) and ( − )/ , respectively, with the same resonant frequency at = 1/√ . Further, according to the transmission spectrum ( Figure 1b) and the characteristic wave impedance ( Figure 1b) that is obtained by a retrieval procedure based on scattering coefficients [51,52], the series LC configuration in Figure 2a fundamentally fits well in terms of the characteristic impedance, i.e., ( − )/ , which has a zero point at . Therefore, a series LC circuit based on the TL circuit A in Figure 2a is chosen to model the resonance of the single-layer metasurface shown in Figure 3. Constructed with symmetrically and periodically arranged elements, such a metasurface shows an ideal geometrical isotropy and polarization independence. So, two basic lumped equivalent circuits with similar band-stop transmission spectra can be used for simple modeling: the parallel and series configurations of LC circuits, as shown in Figure 2, according to the transmission line (TL) theory [51]. meta-atom cell ( Figure 1a). Typically, metasurfaces with such simple capacitive cells of metal elements act as frequency selective devices with transmissive band-stop spectra, as shown in Figure 1b, calculated by the finite element method (FEM). Constructed with symmetrically and periodically arranged elements, such a metasurface shows an ideal geometrical isotropy and polarization independence. So, two basic lumped equivalent circuits with similar band-stop transmission spectra can be used for simple modeling: the parallel and series configurations of LC circuits, as shown in Figure 2, according to the transmission line (TL) theory [51]. Essentially, these two circuits differ in the imaginary parts of their characteristic impedance. Namely, the dominant reactances of the parallel and series LC circuit are − /( − ) and ( − )/ , respectively, with the same resonant frequency at = 1/√ . Further, according to the transmission spectrum ( Figure 1b) and the characteristic wave impedance (Figure 1b) that is obtained by a retrieval procedure based on scattering coefficients [51,52], the series LC configuration in Figure 2a fundamentally fits well in terms of the characteristic impedance, i.e., ( − )/ , which has a zero point at . Therefore, a series LC circuit based on the TL circuit A in Figure 2a is chosen to model the resonance of the single-layer metasurface shown in Figure 3. Essentially, these two circuits differ in the imaginary parts of their characteristic impedance. Namely, the dominant reactances of the parallel and series LC circuit are −jω/ ω 2 − ω 2 0 c and jL ω 2 − ω 2 0 /ω, respectively, with the same resonant frequency at ω 0 = 1/ √ LC. Further, according to the transmission spectrum ( Figure 1b) and the characteristic wave impedance (Figure 1b) that is obtained by a retrieval procedure based on scattering coefficients [51,52], the series LC configuration in Figure 2a fundamentally fits well in terms of the characteristic impedance, i.e., jL ω 2 − ω 2 0 /ω, which has a zero point at ω 0 . Therefore, a series LC circuit based on the TL circuit A in Figure 2a is chosen to model the resonance of the single-layer metasurface shown in Figure 3. meta-atom cell ( Figure 1a). Typically, metasurfaces with such simple capacitive cells of metal elements act as frequency selective devices with transmissive band-stop spectra, as shown in Figure 1b, calculated by the finite element method (FEM).
(a) (b) Constructed with symmetrically and periodically arranged elements, such a metasurface shows an ideal geometrical isotropy and polarization independence. So, two basic lumped equivalent circuits with similar band-stop transmission spectra can be used for simple modeling: the parallel and series configurations of LC circuits, as shown in Figure 2, according to the transmission line (TL) theory [51]. Essentially, these two circuits differ in the imaginary parts of their characteristic impedance. Namely, the dominant reactances of the parallel and series LC circuit are − /( − ) and ( − )/ , respectively, with the same resonant frequency at = 1/√ . Further, according to the transmission spectrum ( Figure 1b) and the characteristic wave impedance (Figure 1b) that is obtained by a retrieval procedure based on scattering coefficients [51,52], the series LC configuration in Figure 2a fundamentally fits well in terms of the characteristic impedance, i.e., ( − )/ , which has a zero point at . Therefore, a series LC circuit based on the TL circuit A in Figure 2a is chosen to model the resonance of the single-layer metasurface shown in Figure 3.  Herein, to take the ohmic loss into account, a resistor R is added into the series with inductance L 1 , and Z o is the vacuum wave impedance. C 1 is the total equivalent capacitance determined by the electric resonance of the single-layer metasurface. It can be modeled as the coupling capacitance C in between the outer rings and inner Jerusalem crosses in parallel with the external coupling capacitance C out between the outer neighboring rings, namely, where C in ≈ πε r ε 0 l 2 / ln(d 1 /t) and C out ≈ πε r ε 0 D/ln[(P − D)/t] are used for good approximations of the geometries [53], and d 1 = D − 2g − 2w 2 − l 1 denotes the internal gap between the cross and square ring, ε r is the relative permittivity of dielectric substrate underneath, t is thickness of the metal ring, and other geometries, l 2 , l 1 , w 2 , D, and g, are shown in Figure 1. Additionally, the inductance of L 1 dominated by outer square rings length can be estimated according to spiral inductance as where µ 0 and ε 0 are the vacuum permeability and permittivity, respectively, of ε r , which is the relative permittivity of the medium in the vicinity, and the other parameters, D, g, l 1 , l 2 , w 2 , and t, are denoted in Figure 1a. As a result, by substituting Equations (1) and (2), the resonant frequency of single-layer metasurface is estimated as [54] Further, the series inductor L s and shunt capacitor C s in Figure 3 denote the dielectric substrate (or Rogers spacer) modeled by transmission line according to the Telegrapher's equations. Specifically, inside the stop band determined by L 1 and C 1 , the inductor of L s and capacitor of C s can be considered as short and open circuits, respectively, due to the relatively small values. Therefore, according to the TL equivalent circuits of two-port networks [54], the characteristic impedance (Z c ), and its refection (R), and the transmission (T) coefficients in the stop band can be analytically approximated as and and respectively, where Z 1 = R 1 + 1/jωC 1 + jωL 1 is the complex impedance of the series LC circuit. Although Equations (1)-(3) tentatively give the equivalent capacitance and the inductance and resonant frequencies, respectively, the lumped elements in Figure 3 can be further adjusted by fitting the resonant frequency, impedance Z c , and the scattering coefficients denoted by Equations (4)-(6) with the counterparts obtained by FEM calculations and impedance retrieval. As a result, spectral responses including the bandwidth and resonant frequency can be accurately modeled and interpreted.

Double Metasurfaces Cascaded with Intercoupling
Further, for a higher degree of freedom in spectral tuning, multiple metasurfaces cascaded in parallel with interlayer couplings are taken into investigations in the following subsection. First, a double-layer setup is constructed by cascading two parallel metasurfaces with one dielectric spacer of Rogers 3003 sandwiched in between, as shown in Figure 4a. To model its spectral behaviors, a lumped equivalent circuit is directly built by cascading two series LC circuits, as shown in Figure 4b.

Double Metasurfaces Cascaded with Intercoupling
Further, for a higher degree of freedom in spectral tuning, multiple metasurfaces cascaded in parallel with interlayer couplings are taken into investigations in the following subsection. First, a double-layer setup is constructed by cascading two parallel metasurfaces with one dielectric spacer of Rogers 3003 sandwiched in between, as shown in Figure 4a. To model its spectral behaviors, a lumped equivalent circuit is directly built by cascading two series LC circuits, as shown in Figure 4b.  Herein, the original series LC circuit is kept to model the top (first) metasurface that faces the incident electromagnetic wave, and the other LC circuit is introduced to model the bottom (second) metasurface. Herein, initial values of capacitances C1, C2 and inductances L1, L2 for both metasurfaces can be estimated individually with Equations (1) and (2). To model their interplay via the dielectric spacer, the same values of Ls and Cs are still used as those in Figure 3b. Further, a mutual inductance of M is introduced between two series LC branches. For a large interlayer gap, when both layers work in independent resonance modes and no mutual coupling exists (i.e., M = 0), the circuit in Figure 4b is simply a parallel combination of two series LC circuits and the spectral behaviors (resonance, transmission, etc.) can be readily obtained in a similar manner as it can in Equations (3)-(6).
However, due to the sub-wavelength gap in our work, strong interlayer coupling exists between two metasurfaces, and a considerable mutual inductance, M, exists between two series LC circuits, as indicated in Figure 4b. In this situation, the total equivalent impedance of two series LC circuits in parallel with mutual inductance can be derived as where = + 1/ + and = + 1/ + denote the initial impedance of either series LC branch. The impedance Z in Equation (7) can be also regarded as the total equivalent impedance of the cascaded metasurfaces. Therefore, with the mutual inductance being considered, the resonant frequencies can now be found when impedance Z in Equation (7) approaches zero. For any fixed value of = , either the positive (additive) or negative (subtractive) mutual inductances, there are constantly two zero points for Z = 0 in Equation (7), which indicate two separate frequencies for the band-stop behaviors of the cascaded metasurfaces.
As a special case, when = = , = = and R1 and R2 are infinitesimal and neglected, respectively, the initial two zero points converge to one. In this situation, either LC series branch has the same equivalent inductance, which is + for the additive mode or − for the subtractive mode of mutual inductance, respectively. Therefore, the equivalent circuit in Figure 4 works at two resonant frequencies of 1/(2 ( + ) ) and 1/(2 ( − ) ), which can be regarded as reversely split from Herein, the original series LC circuit is kept to model the top (first) metasurface that faces the incident electromagnetic wave, and the other LC circuit is introduced to model the bottom (second) metasurface. Herein, initial values of capacitances C 1 , C 2 and inductances L 1 , L 2 for both metasurfaces can be estimated individually with Equations (1) and (2). To model their interplay via the dielectric spacer, the same values of L s and C s are still used as those in Figure 3b. Further, a mutual inductance of M is introduced between two series LC branches. For a large interlayer gap, when both layers work in independent resonance modes and no mutual coupling exists (i.e., M = 0), the circuit in Figure 4b is simply a parallel combination of two series LC circuits and the spectral behaviors (resonance, transmission, etc.) can be readily obtained in a similar manner as it can in Equations (3)-(6).
However, due to the sub-wavelength gap in our work, strong interlayer coupling exists between two metasurfaces, and a considerable mutual inductance, M, exists between two series LC circuits, as indicated in Figure 4b. In this situation, the total equivalent impedance of two series LC circuits in parallel with mutual inductance can be derived as where Z 1 = R 1 + 1/jωC 1 + jωL 1 and Z 2 = R 2 + 1/jωC 2 + jωL 2 denote the initial impedance of either series LC branch. The impedance Z in Equation (7) can be also regarded as the total equivalent impedance of the cascaded metasurfaces. Therefore, with the mutual inductance being considered, the resonant frequencies can now be found when impedance Z in Equation (7) approaches zero. For any fixed value of M = k √ L 1 L 2 , either the positive (additive) or negative (subtractive) mutual inductances, there are constantly two zero points for Z = 0 in Equation (7), which indicate two separate frequencies for the band-stop behaviors of the cascaded metasurfaces.
As a special case, when C 1 = C 2 = C, L 1 = L 2 = L and R 1 and R 2 are infinitesimal and neglected, respectively, the initial two zero points converge to one. In this situation, either LC series branch has the same equivalent inductance, which is L + M for the additive mode or L − M for the subtractive mode of mutual inductance, respectively. Therefore, the equivalent circuit in Figure 4 works at two resonant frequencies of 1/(2π (L + M)C) and 1/(2π (L − M)C), which can be regarded as reversely split from 1/(2π √ LC) due to the different polarities of mutual inductance corresponding to different modes of coupling between two identical metasurfaces. Obviously, the stronger the coupling is, e.g., for a smaller interlayer gap, the larger the mutual inductance M is, and the other two resonant frequencies split farther away from each. In this manner, the stop-band bounded by two resonant frequencies that are tunable by adjusting the interlayer inductance (gap) becomes scalable due to the intercoupling between two metasurfaces. This principle also holds for the general case C 1 = C 2 , L 1 = L 2 for two series LC branches cascaded in parallel with mutual coupling. Other than the mutual inductance M discussed above, the interlayer gap in between also leads to an extra coupling capacitance C m , which can be approximated as [53] C m = 0.2ε 0 ε r1 gD/t 1 (8) where ε r1 and h 1 denote the permittivity and thickness of the dielectric spacer (Rogers 3003) filled in between, respectively. Herein, the interlayer capacitance Cm, which is proven to be much smaller than the intra-layer capacitance is [53], is equivalently incorporated into the single-layer capacitances C 1 and C 2 for simplicity. Therefore, the cascaded model works at distinctly split frequencies due to interlayer coupling, which enables the broadband scalable spectrum. The combined modes of electric resonance and interlayer coupling cause hybridized spectral behaviors of band expanding, as modeled by the LC circuit in Figure 4b. In addition, L s and C s also remain as short and open circuits within the stop band determined by L 1 , C 1 , L 2 , and C 2 . Therefore, given the total equivalent impedance in Equation (7), the spectral behaviors, including the resonant frequencies and the transmission coefficients, can be obtained similarly with Equations (4)-(6) and fitted to the FEM calculation results.

Triple Metasurfaces Cascaded with Intercoupling
Further, for higher degree of freedoms of spectral tuning for broadband regulation, a third metasurface is further introduced for a tri-layer cascaded model, and its equivalent circuit is built in a similar manner. As shown in Figure 5, another series LC circuit is introduced to model the third resonant mode of hybrid spectral behaviors, with intercoupling among each other. Herein, the L 1 and C 1 and L 2 and C 2 branches still model the top and middle layer metasurfaces, respectively, as above. The single-layer equivalent inductance L 3 and capacitance C 3 in the third branch can be estimated in a similar manner using Equations (1) and (2). Note that for large gaps between the adjacent layers of the triple-layer model, the cascaded model of three series LC circuits in Figure 5 still works with no mutual coupling (or inductances). modes of coupling between two identical metasurfaces. Obviously, the stronger the coupling is, e.g., for a smaller interlayer gap, the larger the mutual inductance M is, and the other two resonant frequencies split farther away from each. In this manner, the stop-band bounded by two resonant frequencies that are tunable by adjusting the interlayer inductance (gap) becomes scalable due to the intercoupling between two metasurfaces. This principle also holds for the general case ≠ , ≠ for two series LC branches cascaded in parallel with mutual coupling. Other than the mutual inductance M discussed above, the interlayer gap in between also leads to an extra coupling capacitance Cm, which can be approximated as [53] = 0.2 / (8) where and ℎ denote the permittivity and thickness of the dielectric spacer (Rogers 3003) filled in between, respectively. Herein, the interlayer capacitance Cm, which is proven to be much smaller than the intra-layer capacitance is [53], is equivalently incorporated into the single-layer capacitances and for simplicity. Therefore, the cascaded model works at distinctly split frequencies due to interlayer coupling, which enables the broadband scalable spectrum. The combined modes of electric resonance and interlayer coupling cause hybridized spectral behaviors of band expanding, as modeled by the LC circuit in Figure 4b. In addition, Ls and Cs also remain as short and open circuits within the stop band determined by L1, C1, L2, and C2. Therefore, given the total equivalent impedance in Equation (7), the spectral behaviors, including the resonant frequencies and the transmission coefficients, can be obtained similarly with Equations (4)-(6) and fitted to the FEM calculation results.

Triple Metasurfaces Cascaded with Intercoupling
Further, for higher degree of freedoms of spectral tuning for broadband regulation, a third metasurface is further introduced for a tri-layer cascaded model, and its equivalent circuit is built in a similar manner. As shown in Figure 5, another series LC circuit is introduced to model the third resonant mode of hybrid spectral behaviors, with intercoupling among each other. Herein, the L1 and C1 and L2 and C2 branches still model the top and middle layer metasurfaces, respectively, as above. The single-layer equivalent inductance L3 and capacitance C3 in the third branch can be estimated in a similar manner using Equations (1) and (2). Note that for large gaps between the adjacent layers of the triple-layer model, the cascaded model of three series LC circuits in Figure 5 still works with no mutual coupling (or inductances). However, as it is similar to the double-layer case, strong coupling exists between the adjacent layers when the interlayer gaps reduce down to sub-wavelength level, and mutual inductances are also introduced among three series LC circuits, as shown in Fig-Figure 5. (a) The triple-layers cascaded metasurfaces and (b) the simplified TL equivalent circuit with one series LC circuit shunted plus two parallel LC circuits to mode the hybrid resonance modes with interlayer couplings. However, as it is similar to the double-layer case, strong coupling exists between the adjacent layers when the interlayer gaps reduce down to sub-wavelength level, and mutual inductances are also introduced among three series LC circuits, as shown in Figure 5b. For generality, M 1 and M 2 denote the mutual inductances due to top-middle and middlebottom couplings, respectively. Still, two sets of series inductor Ls 1 , Ls 2 and shunt capacitor C s1 , Cs 2 are inserted in between to model the dielectric spacers inside the gap, as shown in Figure 5b. Similar to the double-layer model, the total equivalent impedance of three series LC circuits in parallel with mutual inductance can be derived as where Z 1 and Z 2 take the same expressions as those in Equation (7), and Z 3 = 1/ 1 R 3 +jωL 3 + jωC 3 denotes the impedance of the third series LC circuit. The resonant frequencies and transmission coefficients are obtainable according to Equations (4)- (6). Obviously, there are three zero points indicated in Equation (9), so three of them are resonant frequencies that constitute and shape the stop-band of the triple-layer cascaded model in this case. These frequencies can be regarded as the resonance modes of each LC branch hybridized with inter-coupling between the adjacent layers. By tuning each square ring length (D 1 , D 2 , or D 3 ) or other geometries individually, the impedances Z 1 , Z 2 , or Z 3 of each LC branch can be regulated in an independent manner according to Equations (1)-(3). Further, by changing the top-middle gap h 1 or the middle-bottom gap h 2 , the mutual couplings or inductances are also adjustable, so are the resonance modes, as well as the spectral behaviors that are determined by the two adjacent layers. As a result, the whole spectrum, including the bandwidth, resonant frequencies, etc., turns out to be scalable in terms of the geometries of the individual layer and inter-coupling effect. Similar to the double-layer case, a smaller (or larger) gap between the adjacent layers for stronger (or weaker) mutual coupling in between cause two adjacent resonant frequencies to split farther (or closer) away from each other.

Results and Discussions
To verify our scheme for the scalable spectral regulation by cascaded metasurfaces with intercoupling, the lumped equivalent circuit models were analyzed and fitted to the numerical calculations by finite element method (FEM). First, the initial values of the lumped elements (capacitance, inductance, and resistance), characteristic impedances, and transmissions (or reflections) of the TL circuit models were calculated using Equations (1)- (9). Thereafter, accurate values were obtained by further fitting the analytical results to those obtained by FEM simulations. Herein, a retrieval procedure based on scattering coefficients was used to extract the characteristic impedance of a single metasurface according to the effective medium theory. Finally, the TL equivalent circuit models are used to instructively tune the spectral behaviors of single and cascaded metasurfaces, e.g., the resonant frequency, bandwidth, etc.

Numerical Verification of the Single-Layer Metasurface
First, for a single-layer metasurface working at the resonant frequency of 50 GHz, all the structural parameters were chosen to be: P = 2.4 mm, D = 1.44 mm, g = 0.12 mm, w 1 = w 2 = 0.12 mm, l 1 = 0.5 mm, and l 2 = 0.28 mm. The metal and dielectric (Rogers 3003) substrate thicknesses are 17.5 µm and 0.762 mm, respectively. After a fitting procedure following the FEM calculation, the lumped elements of the TL equivalent circuit model in Figure 3 were optimized as: L 1 = 1.28 nH, C 1 = 7.8 fF, R 1 = 0.5 Ω, L S = 0.7 nH, C S = 7 fF, and Z 0 = 377 Ω. As shown in Figure 6, the fitted transmission of the TL equivalent circuit model agrees well with the FEM simulation. For spectral regulations through the structural geometries adjustment, the transmission spectra were calculated with varied lengths (D) of rectangular rings (denoted in Figure 1), as shown in Figure 7a. To model and interpret the spectral shifts caused by  For spectral regulations through the structural geometries adjustment, the transmission spectra were calculated with varied lengths (D) of rectangular rings (denoted in Figure 1), as shown in Figure 7a. To model and interpret the spectral shifts caused by variations in the rectangular ring length, the spectra of TL circuit models were also fitted (Figure 7b) with the varied equivalent inductance L 1 , as indicated in Equation (2). Obviously, both the results fit well with the trend of spectral shifts, i.e., a larger ring length leads to a higher inductance and lower resonant frequencies. For spectral regulations through the structural geometries adjustment, the tra mission spectra were calculated with varied lengths (D) of rectangular rings (denoted Figure 1), as shown in Figure 7a. To model and interpret the spectral shifts caused variations in the rectangular ring length, the spectra of TL circuit models were also fit (Figure 7b) with the varied equivalent inductance L1, as indicated in Equation (2). Ob ously, both the results fit well with the trend of spectral shifts, i.e., a larger ring len leads to a higher inductance and lower resonant frequencies. In a similar manner, the other parameters (e.g., l1, l2, and g) also subtly shift resonant frequency of the single-layer metasurface as indicated in Equations (1) and Upon the adjustment of each parameter, the equivalent impedance and transmission efficient of the TL circuit models are obtained using Equations (1)-(6) and fitted to FEM simulated results for accurate predictions. The numerical results of the FEM sim lations and TL models reveal that the outer ring length D dominates over other param ters in shifting the resonance frequency for the single-layer metasurface.

Numerical Verification of Intercoupled Double Metasurfaces
For spectral regulation with more freedom, the cascaded models of m ti-metasurfaces are explored for flexible tuning with more structural parameters, mentioned in Sections 2.  In a similar manner, the other parameters (e.g., l 1 , l 2 , and g) also subtly shift the resonant frequency of the single-layer metasurface as indicated in Equations (1) and (2). Upon the adjustment of each parameter, the equivalent impedance and transmission coefficient of the TL circuit models are obtained using Equations (1)-(6) and fitted to the FEM simulated results for accurate predictions. The numerical results of the FEM simulations and TL models reveal that the outer ring length D dominates over other parameters in shifting the resonance frequency for the single-layer metasurface.

Numerical Verification of Intercoupled Double Metasurfaces
For spectral regulation with more freedom, the cascaded models of multi-metasurfaces are explored for flexible tuning with more structural parameters, as mentioned in Sections 2.2 and 2.3. First, the cascaded model of double-layer metasurfaces was established to work at two split resonant frequencies of around 55 GHz and 45 GHz. The top and bottom metasurfaces were used with the original spacer thickness of 0.762 mm and identical structural parameters of P = 2.0 mm, D = 1.4 mm, g = 0.12 mm, w 1 = w 2 = 0.12 mm, l 1 = 0.48 mm, and l 2 = 0.28 mm. Further, by fitting the ideal transmission spectrum of the TL model to that which was obtained by the FEM calculations (shown in Figure 8), all the lumped elements illustrated in Figure 4 were extracted and optimized as: L 1 = 2.7 nH, C 1 = 4.8 fF, R 1 = 1 Ω, L 2 = 2.7 nH, C 2 = 2.85 fF, R 2 = 1.5 Ω, L S = 1.4 nH, and C S = 3 fF.
In particular, different from the single-layer model, the broadband spectra of double metasurfaces are further scalable by varying the lengths (D 1 or D 2 ) of either of the square rings (top or bottom) or the gap in an independent manner. As predicted in Section 2.2, the cascaded model of the double metasurfaces operates at two split resonant frequencies that can be modeled by two inductively coupled LC circuits, as shown in Figure 8. As shown in Figure 9, the vector electric fields by FEM calculations also confirm the polarity of mutual inductance for two split resonant frequencies. Namely, the lower frequency mode operates with additive mutual inductance, while the higher frequency mode operates with subtractive mutual inductance, which agrees well with our predictions from the cascaded model of two identical metasurfaces. trum of the TL model to that which was obtained by the FEM calculations (shown in Figure 8), all the lumped elements illustrated in Figure 4 were extracted and optimized as: L1 = 2.7 nH, C1 = 4.8 fF, R1 = 1 Ω, L2 = 2.7 nH, C2 = 2.85 fF, R2 = 1.5 Ω, LS = 1.4 nH, and CS = 3 fF. In particular, different from the single-layer model, the broadband spectra of double metasurfaces are further scalable by varying the lengths (D1 or D2) of either of the square rings (top or bottom) or the gap in an independent manner. As predicted in Section 2.2, the cascaded model of the double metasurfaces operates at two split resonant frequencies that can be modeled by two inductively coupled LC circuits, as shown in Figure 8. As shown in Figure 9, the vector electric fields by FEM calculations also confirm the polarity of mutual inductance for two split resonant frequencies. Namely, the lower frequency mode operates with additive mutual inductance, while the higher frequency mode operates with subtractive mutual inductance, which agrees well with our predictions from the cascaded model of two identical metasurfaces.
Furthermore, as shown in Figure 10a, when the top or the bottom ring length increases to be 1600 μm or decreases to be 1200 μm, there is an obvious red or blue shift occurring to the transmission spectra, respectively. Agreeing well with our TL model in Figure 4, any ring length variations distinctly change the equivalent inductance, as well the resonant frequency of that layer or its equivalent LC branch. For the other layer without ring length variations, its frequency also shifts distinctly due to intercoupling or mutual inductance. Specifically, as shown in Figure 10a, when either of the square ring lengths (e.g., D1 for the top ring) decrease from 1400 μm to 1200 μm, the original higher frequency at ~55 GHz (e.g., for D1 = 1400) undergoes a distinct blue shift due to the reduction of its equivalent inductance (L1 for top ring) according to Equation (2). However, for the original lower resonant frequency at ~45 GHz, which is determined by an extra additive mutual inductance (beside the bottom layers L2 and C2), a reduction of either of the ring lengths (D1 or D2) to 1200 μm causes larger deviations between them, and meanwhile, reduced mutual inductance, M. As a result, the lower frequency of 1/(2 ( + ) ) is also increased, which is in good agreement with the discussions in Section 2.2. Furthermore, as shown in Figure 10a, when the top or the bottom ring length increases to be 1600 µm or decreases to be 1200 µm, there is an obvious red or blue shift occurring to the transmission spectra, respectively. Agreeing well with our TL model in Figure 4, any ring length variations distinctly change the equivalent inductance, as well the resonant frequency of that layer or its equivalent LC branch. For the other layer without ring length variations, its frequency also shifts distinctly due to intercoupling or mutual inductance. Specifically, as shown in Figure 10a, when either of the square ring lengths (e.g., D 1 for the top ring) decrease from 1400 µm to 1200 µm, the original higher frequency at 55 GHz (e.g., for D 1 = 1400) undergoes a distinct blue shift due to the reduction of its equivalent inductance (L 1 for top ring) according to Equation (2). However, for the original lower resonant frequency at~45 GHz, which is determined by an extra additive mutual inductance (beside the bottom layers L 2 and C 2 ), a reduction of either of the ring lengths (D 1 or D 2 ) to 1200 µm causes larger deviations between them, and meanwhile, reduced mutual inductance, M. As a result, the lower frequency of 1/(2π (L + M)C) is also increased, which is in good agreement with the discussions in Section 2.2. Apart from the ring length (D1 and D2) tuning, the interlayer gap or dielectric spacer (h1) in between also determines the intercoupling for direct spectral regulations. As shown in Figure 10 (red and blue dashed curves), when the gap shrinks from 762 μm to 506 μm, an obviously blue shift occurs at the original lower frequency, while a red shift occurs at a higher frequency. In good agreement with the predictions of the TL model about mutual inductance discussed in Section 2.2, the larger gap causes two modes to split farther away from each other because the modes with lower and higher frequencies are determined by additive and subtractive mutual inductances, respectively. Namely, the smaller gap causes stronger interlayer coupling and larger deviations in the equivalent inductances of two LC branches in Figure 4. On contrary, as shown in Figure 10 (black curve), when the gap increases to be a very large value (3800 μm), two resonant modes approach each other, and ultimately, merge into one mode because the mutual inductance of two series LC branches or the inter-coupling of two metasurfaces vanishes.
In a similar manner, other structural geometries of l1, l2, and g in both metasurfaces are also tunable for spectral regulations via the double-layer cascaded model. However, we omitted the discussions for such parameters here for simplicity. As a result, by varying the geometries of the top and bottom rings, as well as the interlayer gap, the resonant frequencies and bandwidth are scalable in both a global and local manner.

Numerical Verification of Intercoupled Triple Metasurfaces
For higher freedoms of scalable spectral tuning, the triple-layer cascaded model was verified in a similar manner. As an example, the model was established to work at three split resonant frequencies of 55 GHz, 45 GHz, and 40 GHz, as shown in Figure 11. In this case, the three metasurfaces are customized with slightly different geometries for optimized spectral tuning. For the top metasurface, the parameters of P1 = 2.4 mm, D1 = 1.4 mm, g = 0.14 mm, w1 = 0. 16 (1) and (2), the TL equivalent circuit model is established and perfectly fitted to the FEM results, as shown in Figure 11. Via the fitting procedure, accurate values of the lumped elements denoted in Figure 5 are determined as Apart from the ring length (D 1 and D 2 ) tuning, the interlayer gap or dielectric spacer (h 1 ) in between also determines the intercoupling for direct spectral regulations. As shown in Figure 10 (red and blue dashed curves), when the gap shrinks from 762 µm to 506 µm, an obviously blue shift occurs at the original lower frequency, while a red shift occurs at a higher frequency. In good agreement with the predictions of the TL model about mutual inductance discussed in Section 2.2, the larger gap causes two modes to split farther away from each other because the modes with lower and higher frequencies are determined by additive and subtractive mutual inductances, respectively. Namely, the smaller gap causes stronger interlayer coupling and larger deviations in the equivalent inductances of two LC branches in Figure 4. On contrary, as shown in Figure 10 (black curve), when the gap increases to be a very large value (3800 µm), two resonant modes approach each other, and ultimately, merge into one mode because the mutual inductance of two series LC branches or the inter-coupling of two metasurfaces vanishes.
In a similar manner, other structural geometries of l 1 , l 2 , and g in both metasurfaces are also tunable for spectral regulations via the double-layer cascaded model. However, we omitted the discussions for such parameters here for simplicity. As a result, by varying the geometries of the top and bottom rings, as well as the interlayer gap, the resonant frequencies and bandwidth are scalable in both a global and local manner.

Numerical Verification of Intercoupled Triple Metasurfaces
For higher freedoms of scalable spectral tuning, the triple-layer cascaded model was verified in a similar manner. As an example, the model was established to work at three split resonant frequencies of 55 GHz, 45 GHz, and 40 GHz, as shown in Figure 11. In this case, the three metasurfaces are customized with slightly different geometries for optimized spectral tuning. For the top metasurface, the parameters of P 1 = 2.4 mm, D 1 = 1.4 mm, g = 0.14 mm, w 1 = 0.16 mm, w 2 = 0.1 mm, l 1 = 0.48 mm, and l 2 = 0.32 mm are used, as denoted in Figure 1. For the middle one, the parameters of P 2 = 2.4 mm, D 2 = 1.3 mm, g = 0.16 mm, w 1 = w 2 = 0.16 mm, l 1 = 0.4 mm, and l 2 = 0.28 mm are used. For the bottom metasurface, the parameters of P 3 = 2.4 mm, D 3 = 1.63 mm, g = 0.16 mm, w 1 = w 2 = 0.16 mm, l 1 = 0.56 mm, and l 2 = 0.32 mm are used. By estimating the initial values of all the lumped equivalent elements based on Equations (1) and (2), the TL equivalent circuit model is established and perfectly fitted to the FEM results, as shown in Figure 11. Via the fitting procedure, accurate values of the lumped elements denoted in Figure 5 are determined as L 1 = 3 nH, C 1 = 5.75 fF, R 1 = 2.5 Ω, L 2 = 2.83 nH, C 2 = 3.9 fF, R 2 = 0.9 Ω, L 3 = 3.22 nH, C 3 = 2.15 fF, R 3 = 6.5 Ω, L S1 = 1.2 nH, C S1 = 2.1 fF, L S2 = 1.8 nH, and C S2 = 2.1 fF.  As indicated in Section 2.3, the resonant modes at ~40 GHz, ~45 Hz, and ~55 GHz in Figure 11 are accordingly modeled by the L1, C1 branch, the L2, C2 branch, and the L3, C3 branch together with mutual inductances between the adjacent branches, as shown in Figure 5. As indicated by the tri-layers cascading model, the highest frequency ~55 GHz is determined by the middle layer intercoupling with the top layer; the middle frequency ~45 GHz is determined by top layer intercoupling with the middle layer; the lowest frequency ~40 GHz is determined by bottom layer hybridized with the intercoupling with the middle layer.
Further, as for the scalable spectral tuning by intercoupling, three resonant modes with different frequencies are affected by coupling between the top and middle layers or that between the middle and bottom layers in a different manner. As shown in Figure  12a, at the lowest frequency of ~40 GHz, the FEM-calculated vector electric fields show additive mutual inductances (M2 and M1) for both top-middle coupling (anti-phase) and middle-bottom coupling (anti-phase). At the middle frequency of ~45 GHz, the topmiddle coupling (anti-phase) dominates with an additive mutual inductance (M2) and the middle-bottom coupling (in-phase), which is much weaker than top-middle coupling is, shows a subtractive mutual inductance (M1), as confirmed by the vector electric fields in Figure 12b. At the highest frequency of ~55 GHz, both the top-middle and middle-bottom couplings (both in-phase ones) exhibit subtractive mutual inductances, as shown in Figure 12c. In either case, the middle-bottom coupling is much weaker than the top-middle coupling is, which means the mutual inductance of M2 dominates over M1 (also confirmed by our TL model fitting results). As indicated in Section 2.3, the resonant modes at~40 GHz,~45 Hz, and~55 GHz in Figure 11 are accordingly modeled by the L 1 , C 1 branch, the L 2 , C 2 branch, and the L 3 , C 3 branch together with mutual inductances between the adjacent branches, as shown in Figure 5. As indicated by the tri-layers cascading model, the highest frequency~55 GHz is determined by the middle layer intercoupling with the top layer; the middle frequencỹ 45 GHz is determined by top layer intercoupling with the middle layer; the lowest frequency~40 GHz is determined by bottom layer hybridized with the intercoupling with the middle layer.
Further, as for the scalable spectral tuning by intercoupling, three resonant modes with different frequencies are affected by coupling between the top and middle layers or that between the middle and bottom layers in a different manner. As shown in Figure 12a, at the lowest frequency of~40 GHz, the FEM-calculated vector electric fields show additive mutual inductances (M2 and M1) for both top-middle coupling (anti-phase) and middlebottom coupling (anti-phase). At the middle frequency of~45 GHz, the top-middle coupling (anti-phase) dominates with an additive mutual inductance (M2) and the middlebottom coupling (in-phase), which is much weaker than top-middle coupling is, shows a subtractive mutual inductance (M1), as confirmed by the vector electric fields in Figure 12b. At the highest frequency of~55 GHz, both the top-middle and middle-bottom couplings (both in-phase ones) exhibit subtractive mutual inductances, as shown in Figure 12c. In either case, the middle-bottom coupling is much weaker than the top-middle coupling is, which means the mutual inductance of M2 dominates over M1 (also confirmed by our TL model fitting results). Therefore, other than the individual LC branch, the highest frequency (~55 GHz), middle frequency (~45 GHz), and lowest frequency (~40 GHz) are coupled with the subtractive mutual inductance of M2, the additive mutual inductance of M2 (subtractive mode of M1 is weak and neglected), and the additive inductance of M1, respectively.
First, the top-middle coupling is used for the scalable tuning of the highest and middle frequencies. As confirmed by the FEM results in Figure 13a, reducing (increasing) the dielectric (Rogers 3003) spacer thickness h1 between the top and middle layers causes a distinct red (blue) shift to the middle frequency (~45 GHz) and a blue (red) shift to the largest frequency (~55 GHz), i.e., it drives the two frequencies farther away from (closer to) each other, due to the stronger (weaker) coupling for the larger (smaller) mutual inductance M2. Moreover, the middle-bottom coupling is also used to tune the lowest frequency (~40 GHz) via the additive mutual inductance (M1), since the middle frequency (45 GHz) is mainly determined by the additive mutual inductance of M2. As also confirmed by FEM results in Figure 13b, changing the dielectric gap of h2 for middle-bottom coupling causes distinct red shifts to the lowest frequency (~40 GHz), but with a minor impact on other two frequencies. Therefore, other than the individual LC branch, the highest frequency (~55 GHz), middle frequency (~45 GHz), and lowest frequency (~40 GHz) are coupled with the subtractive mutual inductance of M 2 , the additive mutual inductance of M 2 (subtractive mode of M 1 is weak and neglected), and the additive inductance of M 1 , respectively.
First, the top-middle coupling is used for the scalable tuning of the highest and middle frequencies. As confirmed by the FEM results in Figure 13a, reducing (increasing) the dielectric (Rogers 3003) spacer thickness h 1 between the top and middle layers causes a distinct red (blue) shift to the middle frequency (~45 GHz) and a blue (red) shift to the largest frequency (~55 GHz), i.e., it drives the two frequencies farther away from (closer to) each other, due to the stronger (weaker) coupling for the larger (smaller) mutual inductance M 2 . Moreover, the middle-bottom coupling is also used to tune the lowest frequency (~40 GHz) via the additive mutual inductance (M 1 ), since the middle frequency (45 GHz) is mainly determined by the additive mutual inductance of M 2 . As also confirmed by FEM results in Figure 13b, changing the dielectric gap of h 2 for middle-bottom coupling causes distinct red shifts to the lowest frequency (~40 GHz), but with a minor impact on other two frequencies. Similar as the double-layer model, the spectral behaviors of the triple-layer model are also tunable by changing the square ring length, D1, D2, and D3, of three layers individually. As shown in Figure 14a, increasing (decreasing) the top ring length, D1, from 1400 μm to 1600 μm (1200 μm) causes a distinct red (blue) shift to both the largest fre- Similar as the double-layer model, the spectral behaviors of the triple-layer model are also tunable by changing the square ring length, D 1 , D 2 , and D 3 , of three layers individually. As shown in Figure 14a, increasing (decreasing) the top ring length, D 1 , from 1400 µm to 1600 µm (1200 µm) causes a distinct red (blue) shift to both the largest frequency (~55 GHz) and middle frequency (~45 GHz), but with almost no impact on the lowest frequency (~40 GHz). This agrees well with our TL model above because the top ring length (D 1 ) determines the middle and largest frequencies via the equivalent inductance of L 2 (or L 3 ) and the mutual inductance (M 2 ) between L 2 (top layer) and L 3 (middle layer). Namely, when a larger D 1 induces a larger L 2 according to Equation (2) for a reduced middle frequency, it also leads to weaker intercoupling with the middle layer (due to farther away from D 2 ), and thus, a smaller subtractive mutual inductance of M 2 to L 3 , which ultimately causes a reduction (or red shift) to the highest frequency, and vice versa for the smaller ring length of D 1 .
(a) (b) Figure 13. The FEM calculated transmission spectra of tri-layers model that is scalable by varying (a) the thickness of top-middle coupling dielectric spacer and (b) the thickness of the middlebottom coupling dielectric spacer.
Similar as the double-layer model, the spectral behaviors of the triple-layer model are also tunable by changing the square ring length, D1, D2, and D3, of three layers individually. As shown in Figure 14a, increasing (decreasing) the top ring length, D1, from 1400 μm to 1600 μm (1200 μm) causes a distinct red (blue) shift to both the largest frequency (~55 GHz) and middle frequency (~45 GHz), but with almost no impact on the lowest frequency (~40 GHz). This agrees well with our TL model above because the top ring length (D1) determines the middle and largest frequencies via the equivalent inductance of L2 (or L3) and the mutual inductance (M2) between L2 (top layer) and L3 (middle layer). Namely, when a larger D1 induces a larger L2 according to Equation (2) for a reduced middle frequency, it also leads to weaker intercoupling with the middle layer (due to farther away from D2), and thus, a smaller subtractive mutual inductance of M2 to L3, which ultimately causes a reduction (or red shift) to the highest frequency, and vice versa for the smaller ring length of D1. Further, increasing (decreasing) the middle ring length D2 from 1300 μm to 1500 μm (1100 μm) reduces (increases) the lowest frequency (~40 GHz) distinctly, as well as the other two frequencies, as shown in Figure 14b. Similar to the results in Figure 13a, the middle layer determines the largest and middle frequencies via the top-middle couplings in terms of M2 and the equivalent inductance of L2 (or L3). However, the difference here is the lowest frequency that is also tunable by the middle layer in terms of middlebottom coupling (M1). As an example, when D2 is reduced to 1100 μm, the inductance of L3 decreases to blue shift the largest frequency, and the mutual inductances M2 (topmiddle coupling) and M1 (middle-bottom coupling), both in additive modes as indicated in above discussions, also decrease to cause a blue shift to the middle and lowest frequencies, due to larger differences between D2 and D1 or D2 and D3. Further, increasing (decreasing) the middle ring length D 2 from 1300 µm to 1500 µm (1100 µm) reduces (increases) the lowest frequency (~40 GHz) distinctly, as well as the other two frequencies, as shown in Figure 14b. Similar to the results in Figure 13a, the middle layer determines the largest and middle frequencies via the top-middle couplings in terms of M 2 and the equivalent inductance of L 2 (or L 3 ). However, the difference here is the lowest frequency that is also tunable by the middle layer in terms of middle-bottom coupling (M 1 ). As an example, when D 2 is reduced to 1100 µm, the inductance of L 3 decreases to blue shift the largest frequency, and the mutual inductances M 2 (top-middle coupling) and M 1 (middle-bottom coupling), both in additive modes as indicated in above discussions, also decrease to cause a blue shift to the middle and lowest frequencies, due to larger differences between D 2 and D 1 or D 2 and D 3 .
Finally, for the bottom metasurface, increasing the ring length of D 3 causes a more distinct red shift to the lowest frequency (~40 GHz), with a relatively minor impact on the two higher frequencies, as shown in Figure 14c. As also agrees well with our TL models, since the bottom ring length D 3 mainly involves middle-bottom coupling, which is unable to have a minor impact on the other frequencies. However, when we were reducing the ring length of D 3 to 1430 µm (much closer to D 2 and D 1 ), the middle-bottom coupling increases significantly, and the subtractive mutual inductance of M 1 to causes a distinct red shift to the middle frequency, but with almost no impact on the highest frequency.

Experimental Verifications
As discussed above, by optimally changing the structural geometries of individual layers of cascaded metasurfaces and taking full advantages of interlayer coupling in between them, hybridized resonant modes with split frequencies can be leveraged and tuned in an independent manner for optimal spectral scaling. Therefore, to further experimentally verify the feasibility of the cascaded model for spectral regulations, the optimized devices of both the double-and triple-cascaded metasurfaces were fabricated by using the optimized structural parameters as those obtained in the beginning of Sections 3.2 and 3.3. Herein, Rogers 3003 is embedded as the dielectric spacer in between the cascaded metasurfaces. After sample design and fabrication, the transmission spectra were then characterized in an MM wave experimental testing platform that comprises a pair of horn antennas for MM wave transmitting and receiving and an Agilent vector network analyzer for MM signal analysis. As already shown previously in Figures 8 and 11, the FEM calculations and TL model agree well with each other. Or, in other words, our TL models accurately interpret the spectral behaviors of the cascaded metasurfaces. Therefore, by direct experimental testing of the spectra of two devices that are already designed and optimized according to the TL models, the validities of our TL model for spectral analysis and prediction can be confirmed directly by fitting the measured spectra to the FEM calculated results.
As a result, as shown in Figure 15a,b, the experimentally measured transmission spectra agree well with the FEM calculations (also shown in Figures 8 and 11) for both the double-layered and triple-layered devices. The insets in Figure 15a,b denote the microscopic images of the fabricated device samples, respectively. Obviously, with optimal and scalable spectral tuning, the double-and triple-layer models lead to substantially expanded transmissive stop bands compared with that of the single metasurface in Figure 6. Both the experimental results of double-and triple-layer devices show perfect agreement with the FEM calculations and TL models, which also confirm the validity of our scheme for spectral tuning by simply cascading multiple metasurfaces with similar or the same geometries. red shift to the middle frequency, but with almost no impact on the highest frequency.

Experimental Verifications
As discussed above, by optimally changing the structural geometries of individual layers of cascaded metasurfaces and taking full advantages of interlayer coupling in between them, hybridized resonant modes with split frequencies can be leveraged and tuned in an independent manner for optimal spectral scaling. Therefore, to further experimentally verify the feasibility of the cascaded model for spectral regulations, the optimized devices of both the double-and triple-cascaded metasurfaces were fabricated by using the optimized structural parameters as those obtained in the beginning of Sections 3.2 and 3.3. Herein, Rogers 3003 is embedded as the dielectric spacer in between the cascaded metasurfaces. After sample design and fabrication, the transmission spectra were then characterized in an MM wave experimental testing platform that comprises a pair of horn antennas for MM wave transmitting and receiving and an Agilent vector network analyzer for MM signal analysis. As already shown previously in Figures 8 and 11, the FEM calculations and TL model agree well with each other. Or, in other words, our TL models accurately interpret the spectral behaviors of the cascaded metasurfaces. Therefore, by direct experimental testing of the spectra of two devices that are already designed and optimized according to the TL models, the validities of our TL model for spectral analysis and prediction can be confirmed directly by fitting the measured spectra to the FEM calculated results.
As a result, as shown in Figure 15a,b, the experimentally measured transmission spectra agree well with the FEM calculations (also shown in Figures 8 and 11) for both the double-layered and triple-layered devices. The insets in Figure 15a,b denote the microscopic images of the fabricated device samples, respectively. Obviously, with optimal and scalable spectral tuning, the double-and triple-layer models lead to substantially expanded transmissive stop bands compared with that of the single metasurface in Figure 6. Both the experimental results of double-and triple-layer devices show perfect agreement with the FEM calculations and TL models, which also confirm the validity of our scheme for spectral tuning by simply cascading multiple metasurfaces with similar or the same geometries.

Conclusions
We have demonstrated broadband spectral regulation with flexible scalability by leveraging the intercoupling effects among cascaded metasurfaces in the MMs regime. Doubleand triple-layer cascaded models are constructed by composite meta-atoms composed of Jerusalem resonant crosses and square rings with band-stop behaviors. The resonance modes of periodic composite atoms hybridized with interlayer couplings effects are interpreted and modeled in terms of the simplified TL equivalent circuits to guide the scalable tuning of broad spectral behaviors.
As a proof of concept, the double-layer and triple-layer cascaded devices are constructed and confirmed to scale the narrow stop-band of around 50 GHz for a single-layer metasurface up to a broadened range between 40 GHz and 55 GHz, with an ideal sideband. The experimentally measured results show nicely scaled and enlarged stop-band spectra, which agree well with the predictions of our TL models and FEM calculations. As a result, our methodology based on hybridized resonances with interlayer coupling by cascaded metasurfaces potentially provide an exemplary approach for spectral tuning with broadband scalability in MMs, THz, and related applications.