Excitonic Tunneling in the AB-bilayer Graphene Josephson Junctions

We have considered the AB-stacked bilayer graphene Josephson junction. The bilayers are supposed to be in the charge equilibrium states and at the half-filling in each of the electronic layers of the construction and at each value of the external gate. By considering the interacting bilayers in both sides of the junction and by taking into account both intralayer and interlayer Coulomb interaction effects, we have calculated the normal and excitonic tunnel currents through the junction. The electronic band renormalizations have been taken into account, due to the excitonic pairing effects and condensation in the BLGs. The exact four-band energy dispersions, including the excitonic renormalizations, have been used for the bilayers without any low-energy approximation. We show the degeneracy of the ground state at the zero applied voltage and for different phases of the coherent condensates by showing a dc Josephson current through the junction. The normal and excitonic tunneling currents have been calculated for different gate voltages and for different values of the interaction parameters. The role of the charge neutrality point has been discussed in details.


INTRODUCTION
The existence of the electron-hole bound states for a semimetal with overlapping bands has been postulated long years ago by Keldysh and Kopaev [1], and the prediction about the superfluidity has been given for a condensed excitonic state and ulteriorly, it has been the subject of the intense theoretical studies [2][3][4][5][6][7][8][9][10][11]. Experimentally, the strong evidence of an excitonic insulator (EI) and excitonic Bose-Einstein condensate (BEC) ground states has been shown only in the quantum Hall regime (in a large magnetic field) and under the high pressure in a series of the experimental works on the rare-earth chalcogogenide compounds, transition metal dichalcogenides and tantalum chalcogenides [12]. The exciton condensation was experimentally observed also in quantum Hall bilayers [13,14], in the systems of magnons [15] and cavity exciton polaritons [16,17]. Recently, other solid state systems were proposed as possible candidates for the achievement of the BEC of excitons. It concerns the quantum well heterostructures with the excitons trapped in the cavities of the potential wells [18][19][20][21][22][23][24][25][26][27][28][29] (the structure utilized in these works were doublelayer GaAs/AlGaAs or InAs/GaSb quantum-wells with an electric field applied perpendicularly to the structure), the double layer heterostructures and the bilayers [30][31][32] The excitonic gap formation and the condensation has been examined also in the bilayer graphene structures [33][34][35][36][37][38]. Namely, the bilayer graphene is very promising for the optoelectronic applications due to its unique gate-controllable band structure properties [39]. The imposition of external electrical field can tune the bilayer graphene from the semimetal to the semiconducting state. Nevertheless, the excitonic condensation in the bilayer graphene structures remains controversial in the modern solid state physics because of the complicated nature of the single-particle correlations in these systems [33][34][35][36][37][38]. It has been shown recently [40] that the critical temperature, which describes the transition from the condensate state to the normal state in graphene double layer structure, can be very high due to the extremely small effective mass of excitons. The coherence in exciton BEC condensates survives at the very high temperatures. An analogue conclusion has been drawn in Ref. 38, concerning the bilayer graphene, where the condensate evolution has been analysed as a function of the interlayer Coulomb interaction parameter in the BLG.
Recently, the excitonic condensation has been realized experimentally in the double bilayer graphene heterostructure in the strong quantum Hall regime and by a combination of Coulomb drag and current counterflow measurements [41]. They have also found the evidence of strong interlayer coupling between the graphene layers thanks to the quantized Hall "drag plateau". The zerovalued longitudinal resistance measured there confirms the dissipationless (friction-free) nature of the electronhole condensate state. A quite simple experimental way to observe the excitonic condensate states in the bilayer structures is related to the possibility of engineering of a spatially confined excitonic condensates in the potential traps, and the investigation of the Josephson tunneling effects for excitons [42], related to the tunnelling between two trapped Bose condensates [43] that possess a macroscopic phase coherence. The excitonic Josephson tunnelling effects and thermal transport properties in the electron-hole type double layer graphene junctions, separated by a dielectric layer, have been recently considered in Refs. 44 and 45. In the present paper we study the excitonic tunnelling effects in the tunnel junction based on the AB-stacked bilayer graphene structures. We suppose the presence of the macroscopic phase coherence regime with the well defined condensates phases and amplitudes and we consider only the local on-site interlayer excitonic pairing in each side of the junction. Supposing the electronic bilayers, without the initial optical pumping mechanism, we study the normal and excitonic tunnelling currents through the BLG/I/BLG junction for different values of the external gate voltage applied to the heterostructure. We will assume the half-filling regime in each layer of the BLG structures, even in the presence of the applied gate potential, thus by supposing that not considerable changes of electron density occurs during the adiabatic switching of the external potential. We will show how a finite difference between the phases of the coherent excitonic condensates in the BLG subsystems leads to the excitonic Josephson dc current through the tunnel junction at the zero external voltage. We show also that any finite voltage leads to the ac Josephson current irrespective of the phases of coherent condensates in both sides of the junction. We study the amplitude of the zero voltage Josephson dc current as a function of the interlayer Coulomb interaction parameters in the subsystems. The symmetric and asymmetric interaction cases have been considered straightforwardly. Also, we calculate the normal quasiparticle tunneling current and we show that normal tunneling in the BLG/I/BLG heterostructure is an interaction-protected process and the threshold frequency of the normal tunneling current strongly depends on the values of the Coulomb interaction parameters in the BLGs. We analyse the role of the charge neutrality point (CNP) on the behavior of the normal and excitonic tunneling currents.

A. Description of the Hubbard interactions
We introduce here our model consisting of two Bernal stacking bilayer graphene (BLG) structures separated by a very thick dielectric layer (we suppose that the thickness of the insulating layer is such that we can neglect the quasiparticle scattering and recombination processes in the layer). For the convenience, we will denote by a, b andã,b (and their conjugates a † , b † andã † ,b † ) the annihilation (creation) fermionic operators corresponding to different sublattice sites A, B in the bottom, andÃ,B in the top layer of the left-BLG. Similarly, we denote by c, d andc,d (and their conjugates c † , d † andc † ,d † ) the annihilation (creation) fermionic operators corresponding to different sublattice sites in the bottom, and top layer of the right-BLG. In Fig. 1, we have presented the schematic setup of our BLG/I/BLG junction (here I represents the dielectric layer between the BLGs). We assume here that the layers, which have lattice sites A and B (bottom layer) are biased -V ℓ /2 (with ℓ = L, R), and the layers with the lattice sitesÃ andB (top layers) V ℓ /2, so that the potential difference between the two layers is V ℓ . For the first treatment of such a junction, we suppose the half-filling condition satisfied in both BLG systems, i.e., we suppose that n ℓ ′ ℓ = 1, where the n ℓ ′ ℓ is the total particle number operator in each layer with ℓ ′ = 1, 2 of each BLG with ℓ = L, R. Note, also that we attach the number ℓ ′ = 1 to the bottom layers and ℓ ′ = 2 to the top layers in the heterostructure.
The non-interacting tight-binding part of the total junction-Hamiltonian could be written in the usual form with H ℓ0 , given by where the fermionic operators X ℓ and Y ℓ refer to different sublattice fermions, i.e., for ℓ = L we have X ℓ = a,ã, The intralayer hopping parameter γ 0 is supposed to be the same for both layers and for both BLGs. The fermionic operators P † ℓ and Q † ℓ in Eq.(2) are defined such that for ℓ = L we have P ℓ = b and Q ℓ =ã, whereas for ℓ = R we defined P ℓ = d and Q ℓ =c. The parameter γ 1 denotes the interlayer hopping amplitude between the layers of the BLGs, and it is supposed also the same for both sides of the junction. The variable σ, in all terms in Eq.(2), describes the fermionic spin variable, which takes two values σ =↑, ↓. The last term in Eq.(2) subjects the chemical potential terms with the chemical potentials coupled to the coupled charge density operators in each layer of separated BLGs. Initially, we suppose that the BLG/I/BLG system is in the Grand canonical equilibrium state with the equal chemical potentials in both layers of each BLG structure. Thus, we The charge density operator n ℓ ′ ℓ,σ (r), as it was discussed above, describes the total electron densities with the spin σ in the given layer ℓ ′ = 1, 2 of a given BLG with ℓ = L, R. For example, for the layer ℓ ′ = 1 in the BLG with ℓ = L, we have It is important to notice here that we consider purely electronic graphene layers in both sides of the heterojunction and we do not suppose any initial optical pumping in the BLGs. The fermionic operators in the separate BLG systems satisfy the following anticommutation rules We introduce here also the interaction terms by defining the generalized bilayer Hubbard model in each side of the junction. Namely, we have for the interaction part where H iℓ are the interaction Hamiltonians, for both sides of the construction where n ℓ ′ =1 P ℓ σ and n ℓ ′ =2 Q ℓ σ are the electron densities in the layers 1 and 2 in the BLGs for the P ℓ and Q ℓ -type fermions, introduced above. Then we will use the interaction representation for the fermions [46], in which the time dependence of the fermionic operators is given by the unperturbed Hamiltonian H L0 (V ) + H R0 , where H L0 (V ) is the Hamiltonian of the left side of the junction in the presence of the external gate voltage applied to the BLG/I/BLG. We assume that the gate voltage V (t) drops across the barrier and, in general, the Fermi levels in the left and right-BLG structures will relatively shift by an amount proportional to the potential drop across the junction, i.e.,μ L −μ R ∼ −eV (t). For the reasons that will be clear in the following sections, we have denoted byμ l , l = L, R the exact Fermi levels in different sides of the tunnel junction. We suppose here the half-filling conditions for the total electron densities in each layer of the separate BLG structure and we assume that the influence of the external gate voltage on the charge densities in the layers is infinitesimally small, and the excitonic gap parameter will not get modified by the external perturbation, i.e., δV δn ℓ ′ ℓ δ∆, and consequently δn ℓ ′ ℓ = 0, δ∆ = 0. Such a nonperturbative effect on the electron densities and on the excitonic gap parameter permits to include properly the effect of the applied gate voltage on the excitonic properties in the system, and do not affects the excitonic condensate state in the system. For a more sophisticated case, that evolve the variation of the excitonic condensates states in the junction, one should include the influence of the gate voltage on the excitonic gap parameter and the hal-filling assumption will be failed in this case. Such a treatment is out of scope of the present paper.

B. The tunneling Hamiltonian
The time evolution of the eigenstates is determined by the perturbation term, given via the quasiparticle tunneling Hamiltonian H T . In the following, we will transform the total Hamiltonian of the system H = H 0 + H int + H T into the Fourier space by introducing the Nambu fermionic spinors ψ L k,σ (t) and ψ R k,σ (t) in the left and right sides of the junction. We have The tunnelling matrix Hamiltonian H T will include all possible quasiparticle tunnellings between the BLGs in the junction. Namely, it accounts all sublattice fermions in the hexagonal layers that participate to the tunnelling process. Due to the AB stacking order in the BLGs and the form of the interlayer hopping and interlayer interaction terms in Eqs. (2) and (4), some of the sublattice fermions will not participate to the total excitonic tunnelling (for example the fermions on the lattice sites a,b and c,d), and they will contribute only to the singleparticle normal tunnelling in the junction. Therefore, for a more general case, which includes the tunnellings from all possible fermionic flavours, we will write H T as where t k,p is a matrix element, which describes the transition probability for an electron from the state k (in the left side of the barrier) to a state p (in the right side of it). The total tunnelling current through the BLG/I/BLG junction will be expressed as where Ṅ ℓ=2 R (t) is the expectation value of the rate of change of the electron number operator N ℓ=2 R in the top right layer of the heterojunction, i.e., The expectation value Ṅ R is given by Ṅ R = Tr e −βHṄ R / Tr e −βH , where β = 1/k B T , and H is the total Hamiltonian of the system H. The Heisenberg equation of motion for N ℓ=2 The notation [. . .] − refers to the usual Bose-type commutation rule. A very simple calculation shows that Furthermore, we apply the Ryckayzen [47] transformation for the creation and destruction operators in the left-BLG system. For the top-left and bottom-left layer's electrons, this transformation implies that where X ℓ=2 k,σ (t) =ã k,σ (t),b k,σ (t) in the top layer, and X ℓ=1 k,σ (t) = a k,σ (t), b k,σ (t) in the bottom layer. We can suppose that the phase differences across the junction, between the upper and lower layers in the BLGs, evolve with voltage according to the relations known in the usual Josephson junctions theory [48], and a good reason for this is related to the fact that the total ground state energy of the separate non-interacting graphene layer equals exactly twice of the chemical potential in the layer E ℓ = 2|µ ℓ | with ℓ = L, R [49], which is the case of the usual superconductors according to the Gorkov's derivations of the energy spectrum of superconductors in the microscopic theory [50]. Thus, we have Here, V L and V R are the external gate potentials applied to the left-and right-BLGs, respectively, as it is presented in Fig. 1. The time independent term ∆ϕ 0 /2 in Eq. (13) is given as ∆ϕ 0 = ϕ L − ϕ R and will play the role of the phase detuning parameter in the expression of the excitonic current. Furthermore, ϕ L and ϕ R are the macroscopic phases of the coherent excitonic condensates in both sides of the tunnel junction. Furthermore, we will define the time-dependent phase difference parameter ∆ϕ(t) as where V (t) = (V L (t) − V R (t))/2. Next, we will consider the tunnelling matrix term H T as the small perturbation turned on adiabatically from t = −∞ and which determines the time evolution of the eigenstates in the Heisenberg picture. Therefore, in the first order in H T , we have The infinitesimal, positive constant η = 0 + has been introduced in the integral in Eq. (15) in order to assure the convergence of the integral. Here, the averages over the commutators are referred to the unperturbed Hamiltonian H ′ 0 = H 0 + H int . The time-dependent operators A † (t) and C † (t) in Eq. (15) have been introduced with respect to the reformulation of the tunnelling Hamiltonian in Eq.(7) in terms of the composite tunnel-operators A(t), B(t), C(t) and D(t). Namely, we write C. The normal quasiparticle and excitonic tunneling Next, we calculate the commutators in Eq.(15) and we perform the statistical averaging of four-point correlation functions using the Wick's theorem at the finite temperatures [46]. After the whole averaging procedure, we keep only terms responsible for the normal and excitonic currents in the junction. We obtain for the total tunnelling current the following expression where the first term, in the sum in the right-hand side in Eq. (18), gives the normal single-particle tunnelling current and the second term in Eq. (18) is responsible for the excitonic tunnelling. We obtain Here, we have performed a change of variable of integration t ′ = t − τ and we have introduced the time dependent functions N (t) and F (t) for the normal and excitonic counterparts, which are defined with the help of the Kadanoff-Baym Green's functions [51]. For the normal function N (t), we have and for the excitonic function F (t), we get Th function n ℓ F (ω) = 1/(e βω + 1) (ℓ = L, R) in Eqs. (20) and (21) is the Fermi-Dirac distribution function. We have performed the Fourier transformation of the Kadanoff's Green's functions [51] in order to obtain the expressions of N (t) and F (t), given in Eqs. (20) and (21) and, also, we have introduced the single-particle spectral functions A ℓ X ℓ (k, ω) (X =ã,b for ℓ = L and X =c,d for ℓ = R) and excitonic spectral functions . We have supposed in Eqs. (20) and (21) a simple form of the tunnelling probability t k,p = δ k,p . For a simple treatment, this approximation is sufficient to consider the excitonic effects.
The normal and anomalous (or excitonic) spectral functions in the separate BLGs with the presence of the excitonic pairing interaction have been discussed by us in Ref. 38, where the full four-band theory has been developed without low-energy assumption near the K-point in the Brillouine zone. From the form of the total fermionic action, derived there, it follows that the normal spectral functions in different layers in the left-BLG are interconnected: The same is true also for the right-BLG system. Moreover, the explicit analytical expressions of the normal spectral functions in bottom layer of the left-BLG are [38] A and the similar expressions can be written for the right-BLG structure ℓ = R. The spectral functions in the upper layers in both BLGs could be obtained after the relations in Eq. (22). The k-dependent coefficients α ℓ ik and β ℓ ik are given in Ref. 38. The excitonic energy dispersion parameters κ ℓ ik in Eq. (23), which define the band structures of the interacting BLGs in both sides of the tunnel junction, are given by [38] κ ℓ 1,2k = − As we have mentioned above, the intralayer and interlayer hopping parameters are supposed the same in both sides of the junction. It is particularly important to underline here the role of the bare chemical potentialsμ ℓ appearing in the expressions of the band structure parameters in Eq. (24). Indeed, they are playing the role of the exact Fermi energies in the BLG structures as it was pointed out in Ref. 38. For the anomalous (or excitonic) spectral functions in Eq. (21) we have [38] A (25) and the parameters γ ℓ ik are given in Ref. 38. Furthermore, the normal and excitonic tunnelling currents can be simply expressed analytically after the Werthamer spectral decomposition [52] and by supposing simultaneously the case of the constant gate voltage, i.e., V (t) = V = const. We have and Here, we have used the definition of the phase difference function in Eq. (14). The total tunnelling current through the tunnel junction will be where ω 0 is the field-frequency: ω 0 = e V , and the normal tunnelling is equivalent to the single-particle tunnelling term I n (V, T ) = Im(N (iη + ω 0 )), while the coherent Josephson tunnelling of excitons is given by the last two terms in Eqs. (28). Thus, the total excitonic tunnelling current is and Finally, after some calculations, we get for the normal tunnelling current the following expression Here, we have introduced a product-Θ function Θ T (x, y, z), in the following way where each multiplier, in the total product, is given as a single θ-Heaviside step function. Furthermore, the interaction dependent parameters a 1 , b 1 , c 1 and d 1 , in Eq. (30), are defined as Next, the frequency-dependent parameters ∆ ω0 , ∆ ′ ω0 ,∆ ω0 and∆ ′ ω0 , introduced in Eq.(30), form a detuning matrix∆ and for each component we have where the detuning frequencies ω ′ 0 ,ω ′ 0 , ω ′′ 0 andω ′′ 0 depend explicitly on the difference between the Fermi energies in different sides of the tunnel junction. Namely, we get Next, eight parameters x i with i = 1, ...8, have been introduced in Eq. (30), which are defined as whereas, the function X(x, y, z) is defined as X(x, y, z) = The density of states (DOS) function ρ(x) in Eq.(30) appears after transforming the k-summation in Eq. (20) into the integration over the continuous variable, i.e., k . . . = dxρ(x).... The DOS, in the non-interacting graphene layer, is defined as where γ k is the band dispersion in the non-interacting single graphene sheet, i.e., The parameter d, in Eq. (38), refers to the carbon-carbon distance in the graphene layers. Beyond the Dirac's approximation, the DOS can be analytically expressed [39,53] as where K(x) is the Elliptic integral of the first kind [54] K(x) = π/2 0 dt/ 1 − x 2 sin 2 t . The function Λ(x), in Eq. (39), is given by [53] Λ Furthermore, the functions f (x, y, z) in the denominators in Eq.(30) read as: Concerning the excitonic part of the tunnelling current, the Josephson current I J1 (V, T ), in Eq. (29), is given as: Furthermore, we will calculate numerically the principal value P.V., in Eq.(42), by using the explicit expressions of the excitonic spectral functions, given in Eq. (25). For the Josephson current I J2 (V, T ), we get Again, transforming the summation over the wave vec-tors into the integration, we obtain for the Josephson current I J2 (V, T ): In the following section we will analyse numerically the obtained expressions for the normal and excitonic tunneling currents in the system.

THE NUMERICAL RESULTS AND DISCUSSIONS
A. The normal quasiparticle tunneling In Fig. 2, we have studied the evolution of the normal quasiparticle tunnelling current (given in Eq.(30)) as a function of the applied gate voltage V and for two different limits of the right-BLG interlayer Coulomb interaction parameter W R : the weak and intermediate regime, starting from W R = 0 (solid black curve), W R = γ 0 (solid blue curve), W * R = 1.25γ 0 (bold-dashed darker-green curve) and W R = 1.5γ 0 (dashed darker-yellow curve) (the value W R = 1.5γ 0 is chosen very close to the CNP value W C R = 1.48999γ 0 ) and the other, high interaction limit, when W R = 1.8γ 0 (dashed darker blue curve), W R = 2γ 0 (solid green curve), W R = 3γ 0 (dot-dashed darker red curve) and W R = 5γ 0 (solid red curve). We denoted by W * R the value of the interlayer Coulomb interaction parameter at which the excitonic gap parameter is maximal: The interlayer interaction parameter in the left-BLG is fixed at the value W L = 0.5γ 0 , for all curves in Fig. 2. We see first of all, in Fig. 2, that the normal quasiparticle tunneling through the BLG/I/BLG heterojunction, accompanied with the excitonic pair formations in the system, is a threshold process, and the threshold frequency ω 0 = (e/ )V of the external field depends on the relative values of the Coulomb interaction parameters W L and W R at different sides of the construction. We also see that when augmenting the parameter W R in the interval W R ∈ [0, 1.5γ 0 ], the curves, corresponding to the positive part of the current function I n (V, T ), are shifting into left and the intensity of curves is increased. In turn, the threshold values of the normal tunneling current are also shifting left. Starting from the upper bound (UB) critical CNP value of W R (i.e., W R ≥ W C R (U B)), related to the upper bound solution of the chemical potential in the BLG (see in Ref. 38), the normal tunneling current is shifting right, on the V -axis. Nevertheless, the same is not true for the negative part of the tunnel current. We see, in Fig. 2, that all curves of the negative part of the quasiparticle tunneling current are displacing to right when increasing the interlayer Coulomb interaction parameter in the interval W R ∈ [0, 5γ 0 ]. Only a large jump of the threshold voltages occurs when passing across the lower bound (LB) CNP value W C R (LB), related to the lower bound chemical potential at the CNP. We observe also that for a very large disbalance between the values of the parameters W L and W R , the additional low-  frequency peaks appear in the positive part of the current spectrum, and the threshold frequency values of V are gradually decreasing in these cases. We see also that the amplitudes of the low-frequency peaks are increasing with W R . It is interesting to note that for W R = 5γ 0 (see the solid red curve in Fig. 2), the threshold frequency ω 0 in the positive part of the normal current is of order of ω 0 ∼ 2γ 1 = 0.256γ 0 = 0.76 eV. Another important observation in Fig. 2 is related to the formation of the 4-peak like structures in the spectrum of the normal tunneling current (this is due to the strong excitonic excitations in the 4-band structure of the BLGs). Namely, for small values of W R , the spectrum of the normal current is step-wise, which furthermore transforms in to the 4-peak structure at the intermediate values of W R (W R = 1.2γ 0 , W * R = 1.25γ 0 , W R = 1.5γ 0 and W R = 1.8γ 0 ). This is more apparent in Fig. 3, where we have chosen very close values of W R , in order to demonstrate the gradual formation of the current 4-peak structure. We see in Fig. 3 that for W R = γ 0 the structure of the tunneling current is half-stepwise with the well formed 2 peaks in the excitation spectrum and at the relatively high values of the gate potential. When slightly augmenting W R (see the curves at W R = 1.2γ 0 and W * R = 1.25γ 0 ), the current steps become more pronounced and, at the value W * R = 1.25γ 0 , the resonant tunneling peaks appear at the place of the current steps. Remember that  at W * R = 1.25γ 0 the excitonic gap parameter attains its maximum value (see in Ref. 38). Moreover, for higher values of W R , the 4-peak structure remains present and, additionally, the low-frequency resonant tunneling peaks appear in the electron tunneling spectrum. We relate the 4-peak structure to the high energy strong resonant tunneling of single electrons, and the presence of the very large tunneling threshold is a direct consequence of the EI state in the bilayer graphenes. Thus, in order to do the tunneling, the electrons must break their contribution to the excitonic insulator state. Contrary, the additional low-frequency peaks at the large values of the parameter W R , are related to the excitonic condensate states in the BLGs. This is the manifestation of tunneling coming from the coherent excitonic condensates states in the system. Vis a vis the high interaction values of W R , accompanying the low-frequency peaks, coming from the excitonic condenates states, we can conclude that the excitonic insulator state and the excitonic condensates states in the system are not identical. These are two different states of matter and the EI state does not survives at the high values of W R , considered in Fig. 2. This observation is in complete agreement with the ideas retrieved in Ref. 38. This statement is also in agreement with the recent work in Ref. 44, where it has been shown that the excitonic Josephson current becomes extremely small before the EI state breaks down. A detailed analysis of the role of the CNP point W C = 1.48999γ 0 , for the case of interaction balanced BLGs, is given in Fig. 4. Particularly, in the upper panel in Fig. 4, we consider the equal values of W L and W R , below the upper bound critical value W C L (U B) = W C R (U B) ≡ W C (U B) (the latest was also considered in the figure). We observe that the tunneling spectrum is perfectly symmetric in this case with respect to the origin. When augmenting the interaction parameter up to the upper bound solution at the CNP point, the positive part of the tunneling spectrum is shifted left (see the upper panel in Fig. 4), while for higher values of W R (W R = W L ≥ W C (U B)), the spectrum is shifted to right (this behavior is presented in the middle panel, in Fig. 4). In the lower panel, in Fig. 4, we have presented the smooth passage of the tunneling spectrum when crossing the CNP point. Both, lower bound and upper bound curves of the tunneling current have been considered there. It is worth to mention that the electronic band structure of the BLGs is doubly degenerated at W C because of the chemical potential solutions in the system (see in Ref. 38). We see in the bottom panel in Fig. 4 that the left outermost curve is the UB normal tunneling current in the system. When passing from LB to UB, the tunneling current spectrum still shifted to the left, while for W > W C (U B) the spectrum is transferring to the right. In Fig. 5, we have presented the temperature dependence of the normal tunneling current for a special antisymmetric choice of the Coulomb interaction parameters W L and W R : W L = 0.5γ 0 and W R = γ 0 . As it is clear from the picture, the amplitude of the normal current decreases with increasing the temperature and also leads to the partial suppression of the threshold values of gate voltage, thus promoting a more flexible tunneling of the normal electrons.

B. The excitonic Josephson tunneling
The time dependence of the total excitonic Josephson current through the BLG/I/BLG heterostructure, given in Eq.(29), above, is evaluated numerically in Fig. 6, for different values of the applied gate potential, interaction parameters and condensates phase difference ∆ϕ 0 = ϕ L − ϕ R . When calculating numerically the principal value in Eq.(42), three singular points (−1; 0; 1) of the integrand have been straddled properly. For the left-BLG, we have W L = 0.5γ 0 in all panels in Fig. 6, while for the right-BLG we have chosen three different values W R = 0; 0.5γ 0 and W R = γ 0 , from top to bottom panels in Fig. 6. At V = 0 and ∆ϕ 0 = 0, no tunneling current flows in the system (see the solid black lines on the time axis with the holly triangular plot-markers), while a dc current appears at V = 0 if the phases of coherent excitonic states, in different sides of the junction, differ by ∆ϕ 0 = π/2 or ∆ϕ 0 = −π/2 (see the solid black lines in Fig. 6). Remarkably, the dc excitonic Josephson current changes the sign when changing the sign of ∆ϕ 0 ,  i.e., ∆ϕ 0 → −∆ϕ 0 ⇒ I Exc (t) → −I Exc (t). Such a finite dc current in the system suggests the degeneracy in the ground state of the system, i.e., the U(1) symmetry, and the presence of the excitonic condensates. Moreover, an ac excitonic current appears in the junction for any finite value of V (V = 0.5γ 0 and V = γ 0 , in the picture), even for the case ∆ϕ 0 = 0 (see the blue dotted lines in Fig. 6 with the square plot-markers). The additional phase difference ∆ϕ 0 = π/2 only amplifies the excitonic tunneling current amplitude and leads to a phase shift (see the black dashed lines in Fig. 6). We observe also that the amplitude of the tunnel current decreases when increasing the parameter W R (see the top panel with W R = 0 and the middle panel with W R = 0.5γ 0 , in Fig. 6). In the bottom pannel in Fig. 6 we have chosen larger value of the applied gate potential V = 1.5γ 0 , and W R = γ 0 . We see that the oscillations of the tunnel current are multiplied in this case within the same time interval. This effect is clearly seen in Fig. 7, where different values of the applied gate voltage are considered straightforwardly, and the interaction parameter W R is fixed at the value W R = 0, in correspondence with the plots of the excitonic Josephson current in the upper panel in Fig. 6. We see that when multiplying the value of V by an integer number V ′ = nV , where n = 1, 2, 3, 4, we have for the current wavelength λ ′ = λ/n, thus a relation of type λV = const, between the applied gate potential and the current wavelength, emerges naturally. In Fig. 8, the same function I Exc (t) is presented for the case of the fixed applied gate potential V = 0.5γ 0 , and for different values  of the right-BLG interlayer Coulomb interaction parameter W R , below the critical CNP value W C R = 1.48999γ 0 . The condensates phase difference is fixed at the value ∆ϕ 0 = π/2. We see, particularly, that when augmenting the interaction parameter W R (and keeping at the same time W L fixed at W L = 0.5γ 0 ) the amplitude of the excitonic tunnelling current is decreasing considerably. The zeros of the current function do not shift in their positions, contrary to the case presented in Fig. 7, where the shift of the zeros is caused by the reduction of the current wavelength in the junction. Thus, we realize that the interlayer Coulomb interaction in the right-BLG affects only the current amplitudes, while the changes of the applied gate potential modify principally the frequency of the excitonic current and have not a significant effect on the current amplitudes. Next, in Fig. 9, we have shown the time dependence of the excitonic tunneling current for the values of W R above the UB charge neutrality point W C R (U B). Contrary to the case, given in Fig. 8, the amplitude of the Josephson tunneling current is increasing with W R . This result is related again to the behavior of the chemical potential and Fermi energies in the BLGs (see in Ref.38). In order to compare the results in Figs. 8 and 9, we kept the curve for W R = 0, in both cases. In Fig. 10, we have presented the excitonic Josephson current-voltage dependence on time at T = 0 and for the balanced values of the interaction parameters W L and W R . We put W L = W R = 0.5γ 0 and we have considered the cases when ∆ϕ 0 = 0 or ∆ϕ 0 = 0. We observe, in the upper panel in Fig. 10, that for the case ∆ϕ 0 = π/2  , for V = 0.5γ0 and for different values of the right-BLG interlayer Coulomb interaction parameter WR > W C . The value WR = 0 is also considered as a reference. The parameter WL is fixed at the value WL = 0.5γ0 and the condensates phase difference is fixed at ∆ϕ 0 = π/2. Zero temperature case is considered here. .
we have a finite excitonic Josephson tunneling current through the system at V = 0. At V = 0, we have principally the same qualitative behavior of the tunneling current for both cases ∆ϕ 0 = 0 and ∆ϕ 0 = π/2 apart the situation at t = 0, when there is a large threshold of V at the zero phase difference (see the solid black curve in the lower panel, in Fig. 10). Contrary, for ∆ϕ 0 = π/2 and at t = 0 the excitonic current has been developed in the system (see the solid black curve, in the upper panel, in Fig. 10). The I − V characteristics of the excitonic Josephson current, for the case of the equal values of the interlayer interaction parameters, i.e., W L = W R = W , is shown in Fig. 11. Different values of W are considered in the picture: W = 0.5γ 0 (solid black curve), W = 0.8γ 0 (solid blue curve), W = γ 0 (solid yellow curve), W = W * = 1.25γ 0 (solid red curve), W = 1.3γ 0 (dot-dashed darker green curve), W = 1.4γ 0 (largedashed green curve), and W = W C (LB) = 1.48999γ 0 (bold-dashed darker yellow curve). We observe that up to the value W * = 1.25γ 0 , which corresponds to the maximum of the excitonic gap parameters in the BLGs, the amplitude of the excitonic Josephson current is increasing (including the dc values at V = 0). Furthermore, for W > W * , the amplitudes are continuously decreasing, for W up to the LB CNP value W = W C (LB) = 1.48999γ 0 . The further increase of W , above the LB CNP, leads to a drastic decrease (of about of one order of magnitude) of the excitonic tunnel current amplitude and this is shown in Fig. 12. It is important to concentrate on another principal difference between the results presented in Figs. 11 and 12. This concerns the first deepest minima of the excitonic Josephson current. In Fig. 11, those minima appear for relatively small voltages, and the positions of the positive (negative) minima are shifting to right (left), for W ≤ W * = 1.25γ 0 . Then, with further increase of W in the interval W * < W < W C (LB) = 1.48999γ 0 , they are slightly shifting to left (right). On the other hand, the first minima for the curves in Fig. 12, for W C (U B) ≤ W ≤ 3γ 0 , appear for very large values of V , and the positions of positive (negative) minima are continuously shifting to right (left) in this case. Therefore, the very large values of V are very promising to observe the ac excitonic Josephson current in the system. Experimentally, this could be achieved with the appropriate choice of the left and right gate voltages, which will change the interlayer Coulomb interactions in the BLGs, until the expected effect takes place. In Figs. 13, we have presented the excitonic tunneling current for the asymmetric values of the interaction parameters W L and W R and we see how the increase of W R in the right-BLG, above the symmetric value W R = 0.5γ 0 (the parameter W L is fixed at the value W L = 0.5γ 0 for all values of W R ), leads to the tunneling current transfer into right, on the positive axis of V . It is interesting and straightforward to consider the special case of the non-interacting BLGs and the excitonic Josephson current through the junction in this particular case. This result is shown in Figs. 14. We see that although the zero interaction limit, the ex-

CONCLUDING REMARKS
We have calculated the normal quaiparticle and excitonic Josephson tunneling currents in the BLG/I/BLG heterostructure with the half-filled, AB-stacked bilayer graphene structures at different sides of the junction. By using the S-matrix approximation, we have derived the analytical expressions of both types of tunneling currents: normal quasiparticle and excitonic. The role of charge neutrality point has been discussed in details. Particularly, we have shown how the current spectrum is changing when passing through the CNP point of the interlayer Coulomb interaction parameters in different sides of the junction. It has been shown that the normal quasiparticle tunneling in the junction is a threshold process, independently of the values of the interlayer Coulomb interactions in the BLGs and the field frequency threshold could be modified by changing the interaction parameter in one BLG (for example the right), while keeping it fixed, at the same time, in another BLG in the junction. The formation of the resonant tunneling peaks in the normal tunneling spectrum has been analyzed in details by considering the vicinity of the charge neutrality point in the right-BLG. The very large values of the threshold frequency have been obtained for the case of the interactionsymmetric junction W L = W R , below and above the charge neutrality point. Particularly, from the form of   I −V spectrum, in the case of the interaction-asymmetric junction, it is clear that the excitonic insulator state and the excitonic condensate states are two different states of matter in the BLGs. We have shown that the lowfrequency peaks appear from the excitonic tunneling in the condensates regime in BLGs, and at the large values of the Coulomb interaction parameter in the right-BLG (when the excitonic insulator state breaks down). The temperature dependence of the field frequency threshold values has been derived for a very large interval of temperature: from zero up to very high temperatures.
Furthermore, the time dependence of the excitonic Josephson tunnel current has been studied in details for different limits of interlayer Coulomb interaction parameter in the right-BLG, and the excitonic dc effect has been pointed out for the nonzero phase difference ∆ϕ 0 = 0 and V = 0. Ulteriorly, an ac excitonic Josephson effect has been shown in the junction for non-zero values of the applied gate potential. It has been shown that the phase difference between the excitonic condensates in the BLGs, only amplifies the tunneling current when keeping V constant . It has been shown that a relation of type λV = const is valid in the junction, when changing the applied gate potential by the integer values. The time dependence of the excitonic Josephson tunneling current has been analysed by considering different values of the right-BLG interaction parameter and the role of the CNP point has been revealed out. The I − V characteristics of the excitonic Josephson effect has been found for both interaction-symmetric and interaction-asymmetric cases, and the effects of the charge neutrality point have been discussed in details. It has been shown that the very large values of V are very promising to observe the ac excitonic Josephson current in the system. Finally, a particular case of the non-interacting BLGs junction has been considered separately, and the excitonic Josephson tunneling current has been calculated for this case. As a result, the non zero excitonic current is shown present in this case with the amplitudes comparable to that of the amplitudes of I − V spectrum at the high symmetric values of W in different sides of the junction.
We have presented a self-exhausting theory of the bilayer graphene based heterojunction and the results presented here could represent a veritable framework on which the experimental setup could be made and the results could be compared. The theoretical results in the paper are especially important in the context of the long standing problem about the excitonic condensation and pair formation in the BLGs and double layer graphene systems. This is especially important for theoretical understanding of the nature of the EI state and the coherent excitonic condensate states in such systems: their similarities and differences. The Josephson effect studied in the considered system could be furthermore an important construction in order to build the graphenebased quantum interference devices, which will bring a new idea about the ultra-sensitive magnetometers and voltmeters and will be more sensitive than the usual superconductors-based SQUIDs due to the exceptional mobility of the electrons in graphene. The excitonic Josephson junction is also promising in the context of building a new type of ultrafast electronic circuits blocks, which will form a digital logic unit in the modern ultrafast computers and fast electronics.