where \(n\)  is the index of refraction of graphene in the absence of the applied electric field and \(n_2={6.3e-11}{^2/}\)  is the second-order nonlinear refractive index of graphene [1]}. This formula is based on the assumption that TPA is not accounted for in our system, as we designed our device to operate at THz frequencies, namely in a spectral range where TPA in graphene is negligible [2]}, [3]}. In addition, and equally important, the TPA would not affect the band structure of graphene metasurfaces (the frequency of the photonic bands), as it would merely increase the optical loss of the corresponding photonic states. Moreover, the pump power propagating in an optical mode of a PhC with field distribution, \(\mathbf {E}_p(\mathbf {r})\) , can be expressed as follows [4]}, [5]}:
\(P_p=\frac{v_g}{4a}\int _{V_{cell}}\frac{\partial }{\partial \omega }\left[\omega \epsilon _{m}(\mathbf {r};\omega )\right]\left|\mathbf {E}_p(\mathbf {r})\right|^2d\mathbf {r},\) 
