Designer fermion models in functionalized graphene bilayers
N. A. García-Martínez J. Fernández-Rossier





Departamento de Física Aplicada, Universidad de Alicante, 03690 Spain

International Iberian Nanotechnology Laboratory (INL),
Av. Mestre José Veiga, 4715-330 Braga, Portugal







We propose a method to realize a broad class of tunable fermionic Hamiltonians in graphene bilayer. 
For that matter, we consider graphene bilayer functionalized with sp^3 defects that induce zero energy resonances hosting an individual electron each. The application of an off-plane electric field opens up a gap, so that the zero energy resonance becomes an in-gap bound state whose confinement scales inversely with the gap.
Controlling both the distance among the defects and the applied electric field, we can define fermionic models, even lattices, whose hoppings and Coulomb interactions can be tuned.
We consider in detail the case of triangular and honeycomb artificial lattices and we show how, for a given arrangement of the sp^3 centers, these lattices can undergo an electrically controlled transition between the weak and strong coupling regimes.



Designer fermion models in functionalized graphene bilayers
    J. Fernández-Rossier [On leave from Departamento de Física Aplicada, Universidad de Alicante,  Spain] 
    December 30, 2023
==========================================================================================================





The amazingly rich panoply of electronic phases that occur in solid state matter emerges from the interplay of kinetic energy of the electrons and the Coulomb interaction, both electron-electron and electron-ion. These energy scales are defined by the chemical composition and the crystal structure of each material. Physicists have been looking for strategies to create artificial lattices to confine fermions where the geometry and the energy scales can be tuned independently in order to explore emergent electronic phases escaping the dictatorship of chemistry.
Examples of artificial lattices include optical traps for cold atoms<cit.>, arrays of quantum dots<cit.>, and surface adatoms<cit.>. 


None of these approaches has reached yet the point where a systematic exploration of large scale lattices with non-trivial phases can be controlled in the quantum regime.
Optical traps can create artificial lattices with one fermion per site, but the hopping and Coulomb repulsion energy scales are not much larger than temperature,<cit.> a requirement to observe non-trivial quantum behaviour.
Gate defined quantum dot arrays can definitely avoid that problem<cit.>, but the largest quantum dot array reported so far has only 9 quantum dots<cit.>.
Scanning Tunnelling Microscope (STM) has been used to arrange hundreds of carbon monoxide molecules<cit.> and thousands of chlorine atoms<cit.> on a copper surface, defining artificial lattices for electrons in surface states. However, this approach does not allow further electrical control of the carrier density and the electronic states have finite lifetimes on account of their coupling to the substrate. 

Here we propose functionalized graphene bilayer with an electrically controlled band gap as a flexible and tunable platform to define both fermion and spin model Hamiltonians with a variety of 1D and 2D lattices with a highly tunable energy scales and filling factors.
The building block of the model Hamiltonians are in-gap localized states, produced by sp^3 functionalization of carbon atoms in hollow sites (see Fig <ref> a)). This can be achieved by atomic hydrogen chemisorption<cit.>, but there are many other molecules that produce the same effect<cit.>.
Arrays of these point defects can be created with state of the art STM manipulation,<cit.> but self-assembly of sp^3 molecular centers could be envisioned.
Our calculations show that this platform permits to emulate a number of emblematic lattices (triangular, honeycomb, Kagome, rectangular)

with tunable exchange and Coulomb energies, that implement Hubbard and Heisenberg interactions. 
We show that artificial lattices which realize model Hamiltonians can be created in graphene bilayer which allows the exploration of regimes that lead to the emergence of 2D spin-liquid phases, Mott transitions, quantized anomalous Hall phases, and fractionalized Haldane spin chains and correlated superconductivity. 
 
Different aspects of our proposal have been independently verified experimentally. First, several experiments have demonstrated that an electric field opens up a gap as large as 250meV in graphene bilayer<cit.>. Second, the generation of localized zero modes by sp^3 functionalization has been demonstrated by hydrogen chemisorption in graphene<cit.>, together with the lateral manipulation of the hydrogen atoms, with atomic scale resolution<cit.>. Last, but not least, the observation of correlated electronic phases in twisted graphene bilayer<cit.>, 
that is believed to be associated to the formation of arrays of localized states in a moiré pattern<cit.>, shows that graphene bilayers can indeed host strongly correlated phases.
Here we propose to induce an array of states by chemical functionalization, which would permit to control the location of the bounded states. 
 
We consider Bernal stacked graphene bilayers, shown in Fig. <ref> a). We use the standard single orbital tight binding model<cit.> and consider first neighbour intralayer hopping t=-2.7 and first neighbour interlayer hopping γ=0.4 in accordance with the literature<cit.>.

The effect of sp^3 functionalization is included by removal of a site in the lattice<cit.>: the p_z orbital of the carbon atom forms a strong covalent bond with the functional group and is effectively removed from the p_z lattice<cit.>. We revisit first the problem of a single sp^3 defect in graphene bilayer<cit.>, focusing on the influence of an off-plane electric potential difference, denoted ℰ, which has been considered already in the literature<cit.>.
In the upper layer of Bernal stacked bilayer graphene one of the sublattices is connected to the lower layer while the other sublattice lay on hollow positions. We will consider only functionalization of atoms in hollow sites (see Fig <ref> a)). 

At zero electric field, ℰ=0, pristine graphene bilayer is a zero gap semiconductor, with a degenerate conduction and valence band at the K and K' points of the Brillouin zone<cit.>. The sp^3 functionalization creates a zero mode resonance <cit.>, as expected in any bipartite lattice with a gapless spectrum and with a missing site.<cit.> 
The application of an electric field opens up a gap in the pristine bilayer spectrum and, as we show below, the zero energy resonance becomes a bound state.

In order to study the single particle spectra of graphene bilayer with sp^3 functionalizations we resort to Lanczos diagonalization<cit.> of the single-particle tight-binding model, due to the very large number of carbon sites. For instance, the island used in the calculations of Fig. <ref> has N=131772 sites and a corresponding quantum confinement gap of ∼7 which becomes irrelevant when the gap induced by the electric field is much larger. 
The results of Fig. <ref> c) show the evolution of the in-gap state as well as a few conduction and valence states as a function of the off-plane electric field. It is apparent that as the field is increased, a gap opens and the zero mode resonance stays inside the gap, becoming a bound state with a normalized wave function.



Importantly, the extension of the in-gap bound state induced by the sp^3 functionalization can be electrically tuned. In Fig. <ref> d) we show the confinement length of the in-gap state as a function of the electric field calculated as the minimal radius which contains > 98% of the in-gap wave function.
The extension of the in-gap state is controlled by the size of the energy gap, as shown in the inset of Fig. <ref> d). This result that can be obtained analytically as well<cit.>. This is an important resource in order to control the extension of the in-gap states, and hence the overlap among defects. This will be one of the key ingredients to control the energy bandwidth in artificial lattices and the effective Hubbard repulsion.







We now study the single particle states of artificial arrays of sp^3 functionalization. We have considered two geometries: triangular and honeycomb, although many other structures, such as Kagome and rectangular, are possible by functionalization of the hollow sublattice of the top layer. The calculation is carried out using a supercell

with one (triangular) or two (honeycomb) sp^3 defects. As expected, we obtain as many in-gap bands as sp^3 defects in the unit cell. The corresponding in-gap energy bands for triangular and honeycomb artificial lattices with different electric fields are shown in Fig. <ref>, for two different values of the off-plane electric field.

The in-gap bands can be perfectly fitted to tight-binding models for the triangular and honeycomb lattice with hopping up to third neighbours. These energy bands arise because, in the presence of an electric field, the sp^3 in-gap states are no longer sublattice polarized. Therefore, adjacent in-gap states can hybridize, leading to the formation of in-gap bands. 
The size of the unit cell determines the distance between the point defects and the electric field determines the extension of the bound states.
Thus, both the density of defects and the electric field determine the bandwidth which is in the range of 10 for the cases shown in Fig. <ref>.
 








The results of Figs. <ref> and <ref> show that functionalization of gapped graphene bilayer with ordered arrays of sp^3 defects lead to the formation of in-gap energy bands that can be understood in terms of a tight-binding lattice model, where every functionalization corresponds to a site and a relatively short ranged hopping matrix. 
Importantly, the density of carriers can be tuned independently of the size of the graphene bilayer, using dual gating<cit.>.Therefore, the Fermi energy could be moved along the bandwidth, controlling the filling factor, very much like in the case of the flat bands in the magic angle twisted bilayer<cit.>. 

We now consider the effect of electron-electron interactions. For that matter, we model the Coulomb interaction for the functionalized graphene bilayer in the Hubbard approximation. Our main goal is to provide a fair estimate
of the magnitude of the effective Coulomb repulsion for the electrons occupying the in-gap bands of the artificial lattices. We consider the exact solution for the minimal system where hybridization and Coulomb repulsion compete with each other, namely, a dimer at half filling. The energy scales of this problem are controlled by both the distance between the two sp^3 defects and the gap induced by the electric field.

We first solve the single particle problem for two impurities, using again a bilayer island (with up to N=131772 sites) and two sp^3 centers. The diagonalization yields two in-gap states, with energy ϵ±δ/2, where δ is the hybridization splitting of the states. 




Two states ϕ_1/2, localized states around the defects, can be built out of the two in-gap eigenstates, ψ_±:

    |ϕ_1/2⟩ = 1/√(2)(|ψ_+⟩±|ψ_-⟩)

In the following we refer to ϕ_1 and ϕ_2 as site states, but it has to be kept in mind that they are extended states living in hundreds of carbon atoms. 

We now build the many body Hamiltonian for two electrons in two states, ϕ_1,2. The Hamiltonian can be written as the sum of the singe particle part ℋ_0 and the Coulomb interaction part ℋ_U. The first term reads:

    H_0= t̃∑_σ(c^†_1σc_2σ+c^†_2σc_1σ)
 
where c^†_η,σ, with η={1,2} creates an electron with spin σ in the site state ϕ_η and t̃=δ/2

A minimal model for the Coulomb contribution is obtained using the Hubbard approximation and truncating the Hilbert space to just the two in-gap states. When represented in the basis of the two site-orbitals, the Hubbard term can be expressed as the sum of four types of contributions<cit.>:

    H_ U   =   
    (Ũ_1-J/4) n_1↑ n_1↓ +
    (Ũ_2 -J/4)
     n_2↑ n_2↓ 
     -
       -    J S⃗_1·S⃗_2+  V_ pair+ V_12+  V_21

where n_ησ= c^†_ησc_ησ is the occupation operator in the basis of ϕ_1,2, 
S⃗_η= 1/2∑_σ,σ' c^†_ηστ⃗_σ,σ'c_ησ' are the spin operators associated to these quantum states.

Hubbard-like terms, with energy Ũ_η-J/4 describe the Coulomb energy penalty of double occupation of states ϕ_1 and ϕ_2. The energy scales Ũ_η are given by 

Ũ_η= U ∑_i |ϕ_η(i)|^4


where ∑_i |ϕ_η(i)|^4 is the so-called inverse participation ratio (IPR) and it is a metric of the extension of the state<cit.> and U is the atomic carbon on-site repulsion. Because of the symmetry of site states, we have 𝒰_1≃𝒰_2.
Here we take U=|t|=2.7. The second type of term is a ferromagnetic exchange
with 

J= 2U ∑_i |ϕ_1(i)|^2|ϕ_2(i)|^2


which is a metric of the overlap of the two site states. 
The last two terms are density assisted hopping terms, 

    𝒱_21 = ∑_σ n_1σ(
    t_12 c^†_1σ c_2σ
    +t_12^* c^†_2σ c_1σ)

and the pairing hopping term:

    𝒱_ pair=
    Δ c^†_1↑ c^†_1↓ c_2↑c_2↓
    + Δ^* c^†_2↑ c^†_2↓ c_1↑c_1↓


The Hubbard matrix elements that control these Coulomb energies assisted with hopping processes are: 
t_12= U ∑_i |ϕ_1(i)|^2 ϕ_1(i)^* ϕ_2(i)
and
Δ= U ∑_i (ϕ_1(i)^*)^2 ϕ_2(i)^2.
Both t_12 and Δ could be complex numbers.

The evolution of each of these energy scales as we vary both the applied electric field and the dimer separation is shown in Fig. <ref> a-c). At zero field, the site states are sublattice polarized and the only finite energy scale is the ferromagnetic exchange J and the pair hopping Δ. As the field is ramped up, for both polarities, the site states acquire weight on both sublattices, which permits to have finite t̃, Δ, and t_1,2. At the same time, their extension starts to shrink. The competition between these two effects leads both to the non monotonic behaviour of δ, as well as to the decay of J and Δ. Importantly, for a sufficiently large electric field and defect separation, the dominant energy scales are t̃ and 𝒰, showing that the artificial lattice defines an effective Hubbard model with narrow energy bands. 





We now discuss the properties of the multi-electronic eigenstates of ℋ= ℋ_0+ℋ_U.
We can expand these states in the basis of configurations with a well defined occupation of the single particle states |ϕ_1,2⟩. There are 6 such configurations for 2 electrons in 2 single particle states: four open-shell 
many-body states, where there is one electron per site-state and two closed shell states, where one of the site states is doubly occupied and the other is empty. The eigenstates of ℋ have a well defined total spin S that can take two values, S=0 and S=1. The triplet is formed with open shell configurations.
In contrast, the S=0 states combine both open and closed shell configurations.


For the dimer at half-filling, tuning the electric field can lead to three limiting cases with different types of ground state. For zero electric field, the only non-vanishing energy scales are J and Δ. Since the site state wave functions are real, we have J=2Δ, leading to a S=1 ground state, complying with Lieb's theorem for a bipartite lattice<cit.>. 

At finite electric field both 𝒰 and t̃

rapidly overcome the other energy scales, so that the system effectively becomes a Hubbard dimer. This leads to S=0 ground states. Here, we can have two different limiting situations<cit.>. In the strong coupling limit, 𝒰≫t̃≫ J, the ground state is an open shell singlet with an entangled spin state. In the limit of weak coupling, t̃≫𝒰≫ J, the system behaves like a Hydrogen molecule, with two electrons occupying a bonding molecular orbital. The corresponding many-body energy spectra for the weak and strong coupling cases at finite field are shown in Fig. <ref> d).



We now discuss qualitatively on the different electronic phases that could emerge in the limit of large lattices, on the basis of the predictions for a Hubbard model with first neighbour hoppings, although, other terms, such as beyond-first neighbour hopping, direct exchange and other many-body terms are also present. At half filling, the honeycomb and triangular lattices are expected to undergo a transition from a paramagnetic conducting phase to a Mott insulator that can either have magnetic order <cit.> or be spin liquids <cit.>, as the ratio Ũ/t̃ goes above a critical threshold.

Away from half filling many other non-trivial phases have been predicted. For the triangular lattice at 3/4 filling, a non-coplanar spin canted phase with a Chern number 𝒞=2 has been predicted<cit.> giving rise to a quantized anomalous Hall effect. In the case of the honeycomb, chiral superconductivity has been predicted when the Fermi energy is tuned at the Van Hove singularities<cit.>.
The application of perpendicular magnetic field would lead to the quantized regime at much smaller magnetic fields, on account of the larger unit cell. Finally, superconducting proximity effect, induced by lateral superconducting contacts, might bring an additional knob to induce novel electronic phases, including topological superconductivity.<cit.>
Of course, the sp^3 functionalization approach could be used to create zero and one dimensional structures as well.
 
The experimental realization of this type of artificial lattices will face several technical challenges. So far, the manipulation of hydrogen atoms on graphene has been demonstrated with up to 14 atoms<cit.>. In principle, there is no limit on the number of adatoms that can be manipulated<cit.>. For instance, up to 

80000 chlorine atoms on a Cu(111) surface were manipulated, making use of automatic protocols<cit.>.
Importantly, the number of localized states needed for emergent phases to occur
 is definitely smaller than that: the correlated phases observed in magic angle twisted bilayer graphene occur in samples with 10^4 moiré unit cells.

A second important challenge will be the combination of STM manipulation and dual gating. A possible implementation would be to carry out the atomic manipulation on a graphene bilayer placed on a gated substrate and use the STM as a top gate. This would require to have electrodes on the sample to carry out in-plane transport to probe the electronic phases. 

In summary, we have analysed the potential of functionalized graphene bilayer to produce artificial fermion lattices where the band and Coulomb energies can be tuned by means of gate voltages and geometrical control of the defects. The system proposed here combines the advantage of electric control afforded by quantum dots arrays, and the advantage of reversible structural control given by STM atomic manipulation. The proposed platform could be used to explore a very large variety of electronic phases, including Mott-Hubbard transition, magnetically ordered phases, quantized anomalous insulators, spin liquids and superconducting phases.


Acknowledgements.- 
We acknowledge J. L. Lado for fruitful discussions and technical assistance in the calculations.
We acknowledge grants P2020-PTDC/FIS-NAN/3668/2014, financial support Generalitat Valenciana funding Prometeo2017/139 and MINECO Spain (Grant No. MAT2016-78625-C2). 



apsrev4-1


