Size dependence of excited state properties in zigzag graphene nanoribbons

Density-matrix renormalization group results are reported for the excitonic transitions in zigzag graphene nanoribbons within the framework of the simple Hubbard model. This work shows that a one-dimensional quantum lattice model with open boundary conditions is able to capture the many-body features in the low-energy electronic spectra of these systems alongside the size dependence on the ribbon width. The obtained results agree with the trends of charge and spin excitation gaps available from quantum chemistry predictions.


Introduction
The latest decade has seen a dramatically increasing interest into two-dimensional (2D) graphene for its unique properties which make it a promising material for future electronics [1]. However, one of the main limits preventing the application of 2D monolayer graphene in devices such as field-effect transistors is given by the absence of an OFF state due to its gapless electronic structure. Recently, quasi-one-dimensional (1D) stripes of graphene, also known as graphene nanoribbons (GNRs), have attracted significant interest for the possibility of opening energy gaps in the electronic structure of graphene according to the width and the edge geometry of the chosen ribbon structure [2]. Lately, new strategies have been devised for the synthesis of GNRs by attempting to control their size and edge chirality, ranging from the pioneeristic topdown lithographic methods to bottom-up techniques exploiting aromatic precursors as well as chemical and laser unzipping of carbon nanotubes (CNTs) [3]. Beside armchair-edged graphene nanoribbons (AGNRs), widely investigated for their familybehaviour in energy gaps [4], zigzag graphene nanoribbons (ZGNRs) have attracted the focus of researchers for their spin-related electronic structure, which shows a spatial separation on opposite edges of electrons at the Fermi level. This suggests a spin configuration of these electrons with ferromagnetic coupling along each edge but with antiferromagnetic one between opposite sides [5], which would lead to interesting optical properties [6]. The quasi-1D geometry of GNRs enhances Coulomb interactions between electrons and holes, giving rise to excitonic states, which can deeply affect the optical properties of these nanostructures and consequently the design of new devices. In this work the half-filled Hubbard model Hamiltonian mapped to an open-ended 1D quantum lattice model is used to provide a full many-body treatment of the excited-state properties of ZGNRs going beyond mean-field methods [5] and many-body perturbation theory [6]. Simulation results based on density-matrix renormalization group (DMRG) algorithm [7,13] will be shown in order to address the low-energy excitonic fine structure and the size dependence of related optical excitations in these systems.

Quantum lattice models for graphene nanoribbons
Tight-binding (TB) method has been widely used to understand the basic aspects of electronic structure in graphene-based systems such as graphene nanoribbons both within reciprocal [8] and real space frameworks [9]. In the case of AGNRs, TB method has been able to elucidate the size-dependent family behaviour of their electronic band structure. Such findings were obtained by considering the topological equivalence of AGNRs to brick-type lattice models where the hexagonal symmetry of the honeycomb lattice has been neglected, and in turn to two-leg ladder models with uniform TB hopping parameter [9]. This model was successfully employed for addressing the size-dependence of k = 0 excitonic states in these systems within the framework of the simple Hubbard model [4]. Fig.1 shows the structure of a zigzag GNR with the primitive unit cell highlighted on the honeycomb lattice. The system is composed of N carbon lines joining A-and B-type carbon atoms, thus comprising a total number of 2N sites. From Fig.1 it appears that each edge is composed of atoms entirely of the same type, i.e. A (B) type on the upper (lower) edge, unlike AGNRs for which both types are exposed on each side of the ribbon. It was shown that the 2N × 2N k -dependent Hamiltonian can be mapped to that of a 2Nsite open chain with alternating hopping parameters between atomic-like wavefunctions: t and 2tc (k), where c (k) = cos (ka/2). By numerical diagonalization of the related hamiltonian matrix one can obtain energy bands with peculiar features, such as two partially flat degenerate zeroenergy bands between the Dirac points (k = ±2π/3a) and the edge of the Brillouin zone (BZ) with the corresponding states mainly localized at the ribbon edges. Moreover, c (k) = 0 at the BZ edge, so the bands are highly degenerate, whereas at the Dirac point c (k = ±2π/3a) = 1/2, so the quantum lattice model to be considered becomes equivalent to an open-ended chain with uniform hopping t everywhere [9]. This scheme is particularly convenient for performing DMRG simulations applied to the Hubbard model for investigating correlation effects in ZGNRs, since this algorithm best converges for one dimensional systems with open boundary conditions [7]. Figure 1.
Structure of a zigzag graphene nanoribbon with N carbon lines (thick lines) across its width W joining pairs of A-and B-type carbon atoms. The ribbon is infinite along the horizontal direction (dashed lines). The graphene lattice unit vectors are denoted by a 1,2 and the lattice parameter by a.

Method
In order to compute the low-energy excited state properties of GNRs the simple Hubbard model for π electrons has been considered where c † i,σ and c i,σ are electron creation and annihilation operators on first nearest neighbor sites, is the number of electrons on site i with spin σ, t π and U are the nearest-neighbour hopping parameter and the on-site Coulomb repulsion parameter between two electrons with opposite spins, respectively. The parameter U must be considered an effective parameter taking implicitly into account also longer range Coulomb interactions in the limit of static screening [4]. Since a low-intermediate regime for the Hubbard correlation coupling strength (1 ≤ U/t ≤ 2.2) seems plausible for graphene-based systems [10,11], the DMRG calculations reported here were generally performed for the frequently used setting U/t = 1. The many-particle basis set was defined for half-filled systems with n electrons distributed over n = 2N sites and total spin quantum number S z = 0. Therefore, the size of the basis set and the dimension of the matrix to be diagonalized is D = [n!/(n ↑ !n ↓ !)] 2 , where n ↑ and n ↓ are the numbers of spin-up and spin-down electrons, respectively, with n ↑ = n ↓ = N . Since D is already huge for systems with 16 sites (e.g. 8-ZGNR), a numerical solution of the Hubbard Hamiltonian based on Exact Diagonalization (ED) method would be a computationally prohibitive task. Thus DMRG simulations offer a viable alternative to accomplish the investigation of wider ZGNRs, which can be mimicked by 1D systems with open boundary conditions [7]. Moreover, DMRG allows to address correctly quantum fluctuations beyond mean-field approximation, whose validity has been recent matter of debate in graphene literature [11,12]. The ground state |Ψ GS of the Hubbard Hamiltonian mapped to an open chain has been obtained through a variational optimization search algorithm based on matrix product states (MPS), which is a class of 1D tensor network states underlying the DMRG method [13], as implemented in the Matrix Product Toolkit code [14]. An iterative eigensolver based on the Lanczos algorithm was used. The obtained results for the smallest 4-ZGNR were validated against ED calculations performed with the ALPS libraries [15]. The low-lying excited eigenstates within each S z sector were found by a sequential optimization procedure and repeated orthogonalization with respect to the ground state and (whether available) the already obtained excited states lying above it [13]: first the ground state |Ψ GS is obtained, then the lowest excited state |Ψ 1 that is orthogonal to the ground state i.e. Ψ 1 |Ψ GS = 0; after that one looks for a state |Ψ 2 such that both orthogonalization constraints with respect to |Ψ 1 and |Ψ GS are obeyed, and so on. This differs from traditional DMRG approaches implementing the calculation of excited states through multi-state targetting. The MPS context also allows to calculate efficiently overlaps and expectation values of the optical conductivity operator, hence optical matrix elements between the ground state and the obtained low-lying excited states | Ψ n |v ,⊥ |Ψ 0 | 2 , where v ,⊥ is the velocity (fermionic current) operator for light polarization either parallel or perpendicular to the ZGNR edges [4]. Other quantities relevant to transport properties which can be easily addressed by DMRG-MPS are the charge gap ∆ C and spin gap ∆ S given respectively by where L is the number of sites in the chosen lattice model, N e the number of electrons and S 0 the lowest S tot z component of the total spin such that S z = 0 (S z = 1/2) for even (odd) N e .  Fig.2 the trends of one-photon active excitonic transition energies (left) and related intensities (right) are reported against the ZGNR width expressed by the number of sites in the open chain model. The trend of the transition energies show clearly an inverse dependence on the ZGNR width with only the lowest one-photon active transition for light polarized perpendicular to the ribbon edge |ES11 ⊥ displaying a non-negligible oscillator strength. However, the magnitude of the oscillator strength of this exciton is lower than those found in AGNRs for both polarization directions. Whereas in AGNRs two optical transitions were found for parallel light polarization and one for light polarized perpendicular to the ribbon edge, in this case the situation is reversed. This result agrees with the well known π/2 phase shift of the polarization dependence of interband optical transition in ZGNRs with respect to AGNRs (see Refs. 22, 23 in Ref. [4]). Fig.3 shows the excited state spectrum at the Dirac point of 8-ZGNR obtained with parameters U/t = 1 and t = 2.6 eV, displaying an alternating sequence of singlet and triplet states at low energy. Only one triplet dark (optically inactive) exciton is found here just above the ground state |GS , unlike AGNRs and CNTs, which show a few high-spin states above it.
Remarkably, the two-photon active state denoted as |TPA2 is found to have an oscillator strength value comparable or even slightly higher than the one related to the lowest bright onephoton active exciton |ES11 ⊥ . Such result was not found in AGNRs, where the two-photon intensities were always found to be always lower than the oscillator strength of one-photon active excitons. This finding can be relevant for devising optical limiting applications based on the interesting non-linear optical properties of graphene-based systems [2].     with the complete active space approximation (CAS-CI) [10]. Moreover it was verified that E 0 (N e + 1, S 0 )−E 0 (N e − 1, S 0 ) = U . Fig.5 also confirms their prediction related to the opposite trends of ∆ C and ∆ S with increasing Hubbard correlation parameter U/t, which supports the possibility of developing spin filtering applications based on the combination of both insulating and magnetic behaviour of ZGNRs by tuning electron correlations.

Conclusions
In summary, simulation results were reported for the excited state spectra of zigzag graphene nanoribbons based on density-matrix renormalization group algorithm applied to the halffilled Hubbard model within the recently formulated matrix product state formalism. It was shown that the mapping to an open-ended 1D quantum lattice model captures correctly the polarization-dependent optical features of these systems as well as the size-and correlationdependent trends of charge and spin gaps in agreement with recent results obtained from quantum chemistry configuration interaction methods. Further analysis of the variationally optimized many-body states is required to inquire the spatial profile of excitons and spin polarization, especially at the ribbon edges, with the help of appropriate correlation functions, which will be the subject of future work.