As a typical non-metallic conjugated polymer, two-dimensional graphitic carbon nitride (2D g-C3N4) has received increasing attention in heterogeneous catalysis, photoelectric conversion, sensing, and imaging owing to its appealing electronic band structure, unique chemical structure, and excellent optoelectronic properties [1,2]. Furthermore, its chemical composition and band structure (with a bandgap of around 2.7 eV) are tunable. Because of the luminescent characteristics and biocompatibility of g-C3N4, it has the potential to broaden its application in the field of biosensing [3–6]. Currently, various factors limit the biosensing applications of g-C3N4 based on their photoluminescence and electrochemiluminescence properties [7,8]. It is urgent to develop novel structures that take advantage of the excellent properties of g-C3N4 units for biosensing. Terahertz (THz) metamaterial biosensors combine the local near-field enhancement of metamaterials with the inherent benefits of THz waves in the biological field, showing great potential in the trace detection of biological macromolecules [9–13]. THz metamaterial biosensors research has spread across many biomedical fields, including biological solution detection, microbial detection, and tumor cell screening. The concentration change in biomolecules can be detected by measuring the shift in the frequency range of resonant positions [14–18]. Therefore, the sensitivity of those biosensors needs further improvement. Combining hybridization metamaterials and the tunable energy band features of g-C3N4 results in a novel technique for accomplishing biosensing in the THz range. However, single g-C3N4 has significant limitations, especially a high recombination rate of photogenerated electron–hole pairs, which limits its future application. To suppress the recombination of photogenerated electron–hole pairs, g-C3N4 is recombined with other materials, and the synergistic effect between the two materials can enhance their sensing performance.

Graphene is a 2D, single-atomic layer of carbon material with excellent properties for many potential device applications in the THz range [19–21]. However, graphene has a zero bandgap, which limits its application. The fabrication of graphene–semiconductor heterojunction (HTJ) is one of the approaches suggested for tuning the bandgap [22–24]. Recent studies have taken advantage of the unique properties of different 2D materials by forming an HTJ, which is a vertical stack of dissimilar 2D materials held together by van der Waals (VdW) forces. Gr/MoS2 [25], Gr/h-BN [26], and C3N4-MoS2 [27] VdW HTJ forces have been extensively studied experimentally and theoretically [28,29]. They preserve their intrinsic electronic properties and create outstanding opportunities for engineering applications. In addition, vertical stacking of magnetic metal phases with graphene may be beneficial for promising graphene-based spintronic applications. Most current research on VdW HTJs is based on semiconductor optoelectronic devices and there are few reports on the application of VdW HTJs in biosensing, particularly research on combining with THz metamaterials.

In this study, we propose three THz metasurface biosensors hybridizing three different semiconductor active layers for sensing protein molecules, offering a novel strategy for low-concentration biomolecule detection. The designed metasurface comprises D- and C-shaped metals that produce a dipole resonance of approximately 1.4 THz in response to THz waves with varying protein concentrations. The three different semiconductor active layers are g-C3N4, graphene, and g-C3N4//graphene heterojunction. As a proof of concept, the transmission spectra of the three biosensors were measured with increasing concentrations of casein via the electrostatic doping effect. Some obvious dependences are observed, such as frequency shift, increase in phase difference, and decrease in transmission amplitude with increasing casein concentration. These findings indicate that our biosensors have the ability to identify casein concentrations. Our biosensors can detect casein concentrations as low as 3.54 ng/mL. More importantly, when the casein concentration reaches its maximum (1.56×106ng/mL), the frequency shift and phase difference in our measurements reach their maxima, which are up to 251 GHz and 327°, respectively, for the HTJ biosensors, successfully realizing ultrasensitive, multi-dimensional biosensing. Finally, the Silvaco Atlas package was used to simulate the energy band structure and carrier transport of the three biosensors. The simulation results successfully validated the correctness of our proposed theory model. Our proposed technique could be a viable alternative to future THz biosensing technologies.

Figure 1 illustrates an overview of the fabrication process for the three proposed THz metasurface biosensors. First, a 2 μm thick polyimide (PI) film is spin-coated onto a 300 μm thick quartz substrate (1.5cm×1.5cm) (Fig. 1, ①). Second, aluminum microstructure units are fabricated on the PI film using photolithography technology (Fig. 1, ②). Finally, g-C3N4 is spin-coated and a graphene film is transferred onto the metasurface to produce samples 1 (MS@CN) and 2 (MS@Gr) of the THz metamaterial biosensor, respectively (Fig. 1, ③ and ④). g-C3N4 is spin-coated on the surface of sample 2 to obtain sample 3 (MS@HTJ) (Fig. 1, ⑤). The biosensing feature of the three samples was experimentally characterized via an 8f confocal THz time-domain spectroscopy. Casein molecules were chosen as a probe analyte to verify the performance of the three biosensors (Fig. 1). We prepared suspensions with five different casein concentrations (CCs): C1=3.54ng/mL, C2=56.09ng/mL, C3=852.6ng/mL, C4=1.99×103ng/mL, and C5=1.56×106ng/mL. All the experimental data presented in this work are averaged after three tests.

Figure 2 shows the microstructure design and all geometric parameters of the proposed samples. Figure 2(a) shows a photomicrograph of a unit cell of the microstructure. It consists of D and C shapes with 0.2 μm thick aluminum. The geometric parameters for the D shape are R=50μm and r=40μm and for the C shape are r1=20μm, d=6μm. The periodicity is 130 μm. Figure 2(c) shows the Raman spectra of graphene by 514 nm excitation. The figure shows the G (1578cm−1) and 2D peaks (∼2702cm−1), as well as the full width at half-maximum of the G peak (55.2cm−1). These characteristics indicate that graphene has a high quality and is a monolayer [30]. Figure 2(d) shows the X-ray diffraction (XRD) pattern of g-C3N4 powder. The XRD pattern shows diffraction peaks at 28.96° and no other redundant phases.

The proposed samples are simulated using a frequency-domain solver in CST Microwave Studio. The electric and magnetic boundary conditions were along the x and y directions, respectively, as shown in Fig. 2(b). In the simulation, the permittivity values of the quartz substrate, polyimide, and g-C3N4 are 3.84, 3.4, and 7.86, respectively [10,31]. The Drude model is used to describe the conductivity of monolayer graphene. The electric conductivity of aluminum is 3.56×107S/m. The transmission T(ω) was defined as T(ω)=|Es(ω)/Er(ω)|2, where Es(ω) and Er(ω) are the THz electric field amplitudes of the sample and reference, respectively, after Fourier transformation of the THz pulses [32].

The transmission spectra of the three proposed samples were measured and simulated, as shown in Fig. 3. A sharp resonance dip is exhibited in every transmission spectrum at around 1.4 THz. Further, the experiment results agree well with the simulation results except for a minor difference. These variations are due to unavoidable errors in the fabrication and measurement processes. In addition, the surface current distribution of three samples at 1.41 THz was simulated and shown in the insert of the transmission spectrum in Fig. 3. The surface currents in the D shape of the microstructures clearly follow symmetrical directions. This finding confirms that the resonance at 1.41 THz is caused by dipole resonance. Moreover, owing to the different responses of the three samples to electromagnetic resonance, the surface current intensities of the three samples differ significantly under the same conditions, with the MS@HTJ device exhibiting the strongest inhibitory effect on the surface current.

The experimental transmission spectra of samples with different CCs were measured, as shown in Figs. 4(a)–4(c). To investigate the internal sensing mechanism for the three samples, we simulated the transmission spectra of three samples of g-C3N4 with different conductivities, different Fermi levels of graphene, and different analyte refractive indices [Figs. 4(a1)–4(c1), respectively]. Moreover, the simulation in the transmission spectrum agreed well with the experimental measurements.

Figure 4(a) shows the experimental transmission spectrum of the MS@CN sample with different CCs. The transmission amplitude of the resonance dip slightly decreases under CCs, ranging from 0 ng/mL to 1.56×106ng/mL. In general, protein molecules have a positive charge in the aqueous environment. When the protein molecules are dropped onto the surface of g-C3N4 as an analyte, the protein-induced electrostatic doping effect will cause charge accumulation in the g-C3N4 layer. While its conductivity increases gradually, the loss of the incident terahertz wave increases, resulting in low transmission value. The results reveal that the MS@CN sample can accomplish low CC detection, but its performance as a sensor is poor. In other words, the transmission amplitudes at different CCs did not change much, and there was basically no frequency shift, indicating that g-C3N4 alone responded weakly to changes in protein molecule concentration. To demonstrate sensing mechanism, we simulated the transmission spectra of the MS@CN sample under varied conductivities of g-C3N4, as shown in Fig. 4(a1). When the conductivity is increased from 1 to 125 S/m, the transmission amplitudes of the resonance dip decrease. The simulation and experimental results agree well.

Figure 4(b) shows the experimental transmission spectrum of the MS@Gr sample under different CCs. The results are very interesting, the frequency of the resonance dip shifts over a wide range, and a slight change in the transmission amplitude with the CC increases. When positively charged casein is coated onto the graphene surface, the Fermi level of graphene changes, the conductivity of graphene changes, and the resonance dip shifts via electrostatic doping effect. When CC increased from 0 to 56.09 ng/mL, the resonance dip only occurs with a decrease in transmission amplitude. When the CC reached 852.6 ng/mL, the resonance dip shifted over a wide frequency range, and the CC caused a significant change in graphene conductivity. As the CC continued to increase, the resonance dip shifted until the CC reached 1.56×106ng/mL; the saturation concentration and resonance dip shift reached their maximum values. To analyze the internal sensing mechanism, we simulated the transmission spectrum of the MS@Gr sample under varied Fermi levels of graphene, as shown in Fig. 4(b1). When the Fermi level of graphene changed from 0.1 to 0.2 eV, the frequency of the resonance dip shifted. The simulation and experimental results agree well.

Based on the research results of the above two biosensors, we designed and used the above two active materials to form a VdW HTJ and then studied the sensing performance of the MS@HTJ sample combined with the HTJ and the metasurface, as shown in Fig. 4(c). The results were as we expected; with the increase of the CC, the resonance dip shifted to the low-frequency direction, and when the concentration reached 1.56×106ng/mL, the resonance dip stopped shifting. Because protein addition changes the charge carriers concentration at the HTJ surface and thus the conductivity via electrostatic doping effect, the frequency of the resonance dip shifts as the HTJ conductivity changes. For the MS@HTJ sample, we simulated the transmission spectrum with different refractive indices, as shown in Fig. 4(c1). The results show that as the refractive index increases, the resonance dip frequency gradually shifts to a low frequency. The simulation results are consistent with the experimental results. The correctness of our proposed theoretical model has been verified.

To further investigate the sensing mechanism of the three samples, we first explain the mechanism of the signal change in the MS@CN sample using the classical semiconductor conductivity formula σ=|e|(μene+μhnh), where e is the electron charge; μe, μh are the electron and hole mobilities, respectively; and ne, nh are the electron and hole concentrations, respectively. According to the formula, when the positively charged protein is doped to the g-C3N4 surface, it causes the rapid separation and migration of electrons at the g-C3N4 interface, causing an instantaneous increase in the surface conductivity of g-C3N4. It enhances destructive interference of THz waves, ultimately causing a decrease in the transmission values at resonance dip as the CC increases. Figure 5(a) shows that the band structure of g-C3N4 changes under different CCs. Region (i) in Figure 5(a) shows without protein, g-C3N4 is in the inactive state of an n-type semiconductor and the conduction and valence bands are not bent. As the CC increases, it is equivalent to incorporating a positive charge on the surface, while a large number of electrons gather on the surface, causing the conduction and valence bands to bend, as shown in regions (ii) and (iii) in Fig. 5(a).

For the sensing mechanism of MS@Gr and MS@HTJ samples, we explain the frequency shift using the modified perturbation theory [13,17]. The electric field can be represented by two equations-based E(r) to H(r) in electromagnetic dynamics, as follows: ∇×∇×E(r)=(ωc)ε(r)E(r).(1)

By applying the perturbation procedure to Eq. (1), we obtain the following formula for the change in the dielectric environment Δε and angular frequency shift Δω values: Δω=−ω2∫|E(r)|2Δεdr∫|E(r)|2ε(r)dr+o(Δε2).(2)

Consider the case where the electric field decays exponentially along the direction normal to the material; Δω calculated from Eq. (2) scales with the analytes as follows [17]: Δωω∝Δε.(3)

In the case of local confinement, the surface conductivity σg of graphene is represented, including the interband term σinter and the intraband term σintra approximation [33,34]: σg=σinter+σintra=ie24πℏln(2EF−(ω+i/τ)ℏ2EF+(ω+i/τ)ℏ)+ie2kBTπℏ2(ω+i/τ)(EFkBT+2ln(eEFkBT+1)).(4)

Here, ω is the angular frequency, e is the electron charge, kB represents the Boltzmann constant, EF is the Fermi level of graphene, ℏ is the reduced Planck constant, τ refers to the relaxation time of graphene, and T is the temperature in Kelvin. The permittivity of graphene is calculated by [35] εg=1+σiε0ωd,(5)where ε0 is the vacuum permittivity and d is the thickness of graphene (1 nm).

For the MS@Gr and MS@HTJ samples, the surface conductivity σg can be changed when doping the samples with charged protein molecules by changing the Fermi level of graphene in Eq. (4). The dielectric constant εg of the graphene layer has been changed in Eq. (5). According to Eq. (2), because εg has been changed, resulting in a dielectric environment in the altered metamaterial, an increase in CC will increase the frequency shift Δf. In this work, the graphene is synthesized using chemical vapor deposition and its Fermi level is downshifted through chemical doping, forming p-type doped graphene. The initial quasi-Fermi level (EF0) is in the valence band near the Dirac point [36–38]; see region (i) in Fig. 5(b). As the CC increases, the quasi-Fermi level upshifts from EF0 to EF1 [as shown in region (ii) in Fig. 5(b)] and then to EF2 [as shown in region (iii) in Fig. 5(b)], while graphene conductivity gradually decreases. The respective band structures of the MS@HTJ samples without containing casein are shown in region (i) in Fig. 5(c); the HTJ is formed by n-type g-C3N4 and p-type graphene. The band structure at the interface in the MS@HTJ sample is compared to that of g-C3N4 in the curved MS@CN sample. With the increase of CC, the bending degree of the band structure of g-C3N4 is strengthened, while the Fermi level of graphene decreases. When the concentration increases to C4, the energy band of g-C3N4 bends to the point where electrons transfer to the graphene interface, and thereby the Fermi level of graphene upshifts into the conduction band, as shown in region (iii) in Fig. 5(c).

In summary, the three samples demonstrated ultrasensitive, multi-dimensional protein sensing functions. To analyze and compare their sensing performances, we quantified the transmission amplitude, frequency shift, and phase difference at the resonant frequency of each sample, as shown in Figs. 6–8.

The color curves in Figs. 6(a) and 6(b) show the frequency and transmission values of resonance in the three samples at various CCs. For the MS@CN sample, there was little shift in frequency with increasing CC; however, for the MS@Gr and MS@HTJ samples, there was a larger shift in frequency. Using the sample without protein as the initial value, the frequency shift of the resonances Δf=|fi−f0| at different CCs was calculated for the three samples, as shown in Fig. 6(c). As the CC increases from C1 to C5, the maximum Δf of the MS@CN sample is 10 GHz. When the CC reaches C5, the maximum Δf is 150 GHz for the MS@Gr sample. Simultaneously, for the MS@HTJ sample, when the CC increases, the Δf changes from 32 GHz to the maximum value of 251 GHz. Compared with the other two samples, the presence of the HTJ causes the response of the resonance dip and the Δf to demonstrate the largest values. In addition, the transmission value changes are different for the three samples. As shown in Fig. 6(b), with the increase of CC, the transmission value of the MS@CN sample decreases all the time, while the transmission values of the MS@Gr and MS@HTJ samples first decrease and then increase. Using the sample without protein as the initial value, the transmission difference for the resonances T=|Ti−T0| at different CCs was calculated for three samples, as shown in Fig. 6(d). As the CC increased, the ΔT change of the MS@CN sample gradually increased from 0.02 to 0.08, while the ΔT change of the MS@Gr sample first increased to 0.11 and then decreased to 0.09, and the ΔT change of the MS@HTJ sample was the least obvious; the largest ΔT value was 0.05. The MS@CN sample has the most noticeable change in the transmission value, while the MS@HTJ sample has the most noticeable frequency shift.

To demonstrate that the advantages of HTJ for the sensor come from the large-scale carrier transport caused by its internal electric potential, the Silvaco Atlas package was used to simulate the energy band structure and carrier transport of the three samples, as shown in Fig. 7. Silvaco Atlas is a 2D device simulation package used for structural modeling by solving a set of fundamental equations comprising the Poisson’s equation, the continuity equations, the transport equations, and the Fermi–Dirac statistics [39,40]. The electrical parameters of materials used in the Silvaco TCAD simulation are listed in Table 1. The device structures are shown in Fig. 1. In experiment, the protein solution has different positive fixed charges concentrations due to dilution. The interface charge density is specified in our simulation study. Figures 7(a)–7(c) show the energy bands and activation energies (Ea) of g-C3N4 without/with positive fixed charges at the surface. Compared with the flat band condition under Qf=0, the conduction band edge near the front surface is close to the Fermi level (EF) of the g-C3N4 substrate. Compared with surfaces without containing fixed charges, the electron concentration near the surface is much higher than the electron concentration inside [as shown in Fig. 7(c)]. Figure 7(b) indicates that the Ea decreases from 0.357 eV for Qf=+104cm−2 to 0.0398 eV for Qf=+1012cm−2. Figures 7(d)–7(f) show the energy bands and activation energies of graphene without/with positive fixed charges at the surface. Figure 7(d) shows that the work function of graphene is not modulated by positive Qf. Even though an external electrical field is present, the induced electron accumulation and hole repulsion can be ignored [as shown in Fig. 7(f)]. The modulation of the work function of graphene usually depends on its interaction with its counterpart material. Figure 7(e) shows that the Ea of graphene of 0 eV is not affected by Qf. The absence of a bandgap in graphene makes it semi-metallic. Figures 7(g)–7(i) show the energy bands of the g-C3N4-graphene HTJ without/with positive fixed charges at the surface of g-C3N4. For this structure, which is also called the g-C3N4-graphene Schottky junction, the conduction band edge (Ec) near the g-C3N4-graphene interface is close to the Fermi level (EF) of the g-C3N4 substrate, indicating electron accumulation at the interface. If a protein with positive fixed charges is added on the surface of g-C3N4, the Ec bends toward the EF up to a magnitude depending on the density of the positive fixed charges (Qf). The Ea variation tendency as a function of Qf in Fig. 7(h) indicates a high Ec of 0.357 eV with a very low positive Qf, and low Ec values of 0.261 eV, 0.142 eV, and 0.02 eV for Qf=+108, +1010, and +1012cm−2, respectively. The decreased Ea and Ec indicate that the concentration of electrons increases in g-C3N4, resulting in the upshifting of the EF of its adjacent graphene. The above simulation results verify the correctness of our proposed theoretical model. It was the first time that Silvaco TCAD was used to analyze the internal sensing mechanism of THz sensors, which provides a new strategy for the application of HTJs in THz technology.Table 1.

Electrical Parameters of Silvaco TCAD Simulation in This Work

The phase differences (ΔP) of the three samples with and without proteins were obtained for CCs ranging from 3.54 to 1.56×106ng/mL, as shown in Figs. 8(a), 8(b), and 8(c), respectively. It is worth noting that the dependence of ΔP on frequency is quasi-linear. The higher the frequency and concentration, the larger the phase difference, as shown in Figs. 8(d)–8(f) of the 3D phase diagram. These results reveal that the phase is also altered after the introduction of protein, mainly because of the change in conductivity of g-C3N4, graphene, and HTJ. Therefore, ΔP can be considered a biosensing index for the three samples. Furthermore, the MS@HTJ sample has the largest ΔP of the three, measuring 327° at 1.56×106ng/mL. As the CC increases, the ΔP shows an obvious change.

For the MS@Gr sample, ΔP does not change significantly as the CC increases; the maximum achievable ΔP value is 159°. Although the ΔP of the MS@CN sample increases as the CC increases, the overall ΔP is small and the maximum value is 116°. When ΔP is used as a dimension of sensing, the sensing effect of the HTJ is the best. Therefore, our claim that the phase change of the MS@HTJ sample is mainly caused by the HTJ has been verified.

The results show that the proposed three samples can function as a multi-dimensional THz biosensor capable of detecting ΔT, Δf, and ΔP in the presence of proteins.

In this study, three novel terahertz biosensors are designed, fabricated, characterized, and comparatively analyzed. The three biosensors include g-C3N4, graphene, and an HTJ integrated with a dipole-resonant metasurface. As a result, the proposed THz biosensors are demonstrated for ultrasensitive, multi-dimensional sensing of casein molecules with a limit of detection of 3.54 ng/mL by identifying the changes in frequency, phase, and amplitude. When the casein concentration reaches its highest value (1.56×106ng/mL), the frequency shift and phase differences reach their maximum values, which are up to 251 GHz and 327° for the HTJ biosensors. We presented simulation calculations and a theoretical analysis of changes in the conductivity of the active layer via electrostatic doping effect to explain the internal mechanism of the biosensor. As casein concentration increases and the covalent bonding of graphitic carbon nitride, graphene, and HTJ strengthens, the energy band structure changes and the conductivity of the active layer changes. More importantly, the MS@HTJ sample showed the greatest response to the protein. To prove the advantages of the HTJ for the sensor, the Silvaco Atlas package was used to model the three samples’ energy band structure and carrier transport. The simulation results successfully verified the correctness of our proposed mechanism theory model. These findings are promising for designing ultrasensitive HTJ THz meta-sensors.

Acknowledgment. We thank OSA Language Editing Group for editing the English text of a draft of this paper.