## 2 Electrical Response of a Bi-Layer Graphene Nanoribbon THz Detector
### Transport Model along Graphene Nanoribbons, with the Tunneling Effect  
In this paragraph, we extend the above model to include the possibility that the two GNRs may interact among them via the tunneling effect. In this case, we modify the generalized [PERSON]'s law (Equation (2)) taking into account the possible mutual interaction between the two GNRs:  
\[\begin{cases}\hat{J}_{z1}=\hat{\sigma}_{11}(\beta,\omega)\hat{E}_{z1}+\hat{ \sigma}_{12}(\beta,\omega)\hat{E}_{z2},\\ \hat{J}_{z2}=\hat{\sigma}_{21}(\beta,\omega)\hat{E}_{z1}+\hat{\sigma}_{22}( \beta,\omega)\hat{E}_{z2},\end{cases} \tag{17}\]  
where \(\hat{\sigma}_{11}=\hat{\sigma}_{22}=\hat{\sigma}_{s}\) are the self-conductivities given as in Equation (2), whereas the mutual conductivities \(\hat{\sigma}_{12}=\hat{\sigma}_{21}=\hat{\sigma}_{m}\) may be expressed as follows:  
\[\hat{\sigma}_{m}(\beta,\omega)\equiv-i\frac{2e^{2}}{\pi\hbar\
u_{F}X}\frac{ \omega_{i}\omega^{\prime}}{\beta^{2}}\frac{1}{1-\left(\beta\
u_{F}/\omega^{ \prime}\right)^{2}}\hat{\gamma}(\beta,\omega) \tag{18}\]  
where the complex frequency is defined as \(\omega^{\prime}=\omega+i\
u\), and  
\[\hat{\gamma}(\beta,\omega)=\sum_{k=1,2}\frac{\omega^{\prime}}{\beta^{2}}\frac {2\omega_{i}+\left(-1\right)^{k}\left[1+\left(\beta\
u_{F}/\omega^{\prime} \right)^{2}\right]}{\left(2\omega_{i}+\left(-1\right)^{k}\omega^{\prime} \right)-\beta^{2}\
u_{F}^{2}} \tag{19}\]  
In the above relations, a fundamental role is played by the tunneling frequency \(\omega_{t}\), related to the binding energy \(\hbar\omega_{t}\) due to delocalization of the \(\pi\)-electrons for the tunnel effect.  
As a consequence of new constitutive relation Equation (17), the MTL model derived by coupling Equation (17) with [PERSON]'s equation would be expressed through the p.u.l. matrix parametersobtained by adding to the inductance and capacitance matrices the operators \(\mathbf{L}_{nm}(\mathbf{I})\) and \(\mathbf{C}_{tun}(\mathbf{V})\), acting as follows:  
\[\mathbf{L}_{nm}\left\{\mathbf{I}(z,\omega)\right\} = \frac{1}{4\pi}\frac{\omega_{i}\omega^{\ast}}{\omega_{V_{F}}}\Bigg{(} \frac{c}{v_{{}_{F}}}\Bigg{)}^{2}\Theta\mathbf{L}_{M}\int\limits_{\rightarrow}^ {\rightarrow\infty}K_{i}(z-z^{\ast})\mathbf{I}(z^{\ast},\omega)dz^{\ast} \tag{20}\] \[\mathbf{C}_{nm}^{-1}\left\{\mathbf{V}(z,\omega)\right\} = -\frac{\omega_{i}\omega^{\ast}}{2v_{{}_{F}}}\Theta C_{E}^{-1}\int\limits_{ \rightarrow}^{\rightarrow\infty}K_{z}(z-z^{\ast})\mathbf{V}(z^{\ast},\omega)dz ^{\ast} \tag{21}\]  
where \(\Theta\) is the matrix:  
\[\Theta=\begin{bmatrix}1&-1\\ -1&1\end{bmatrix} \tag{22}\]  
The kernels in Equation (20) are defined as  
\[K_{1}(z;\omega^{\ast})=\sin\Bigg{(}\frac{2\omega_{i}}{v_{{}_{F}}}\big{|}z \big{|}\Bigg{)}\exp\Bigg{(}-i\frac{\omega^{\ast}}{v_{{}_{F}}}\big{|}z\Bigg{|} \Bigg{)} \tag{23}\]  
\[K_{2}(z;\omega^{\ast})=\frac{1}{2i}\exp\Bigg{(}-i\frac{\omega^{\ast}}{v_{{}_{F }}}\big{|}z\big{|}\Bigg{)}\Bigg{[}\frac{1}{2\omega_{i}-\omega^{\ast}}\exp \Bigg{(}i\frac{2\omega_{i}}{v_{{}_{F}}}\big{|}z\big{|}\Bigg{)}+\frac{1}{2 \omega_{i}+\omega^{\ast}}\exp\Bigg{(}-i\frac{2\omega_{i}}{v_{{}_{F}}}\big{|}z \big{|}\Bigg{)}\Bigg{]} \tag{24}\]  
Note that the operator \(\mathbf{L}_{tun}\) also adds to the kinetic inductance \(L_{k}\mathbf{U}\) in the modified expression of the resistance matrix \(\mathbf{R}\). Operator's Equations (20) and (21) introduce a spatial and frequency dispersion, which affects the relations between the voltages and the p.u.l. charge and the currents, and the p.u.l. magnetic flux.