Bias induced ferromagnetism and half-metallicity in graphene nano-ribbons

Towards spin selective electronics made of three coordinated carbon atoms, here we computationally propose robust and reversibly bias driven evolution of pristine undoped graphene nano-ribbons(GNR) into ferromagnetic-semiconductor, metal or a half metal, irrespective of their edge configurations. The evolution is a result of a rare ferromagnetic(FM) order emerging among nearest neighbouring(n-n) sites, in positively biased regions in their in-homogeneous bias unit-cells, in attempt to cooperatively minimise on-site Coulomb repulsion and kinetic energy, while maximising localization of electrons at the positively biased sites. The phenomenon appears to be a general property of in-homogeneously biased Coulomb correlated bipartite systems. Consequences are particularly rich in zigzag edged graphene nano-ribbons(ZGNR) due to the contest of bias driven n-n FM order and the inter-edge antiferromagnetic order inherent to ZGNRs, leading to systematic closing of gap for one of the spins, amounting to bias controlled unmissable opening of window for FM-semiconducting and half-metallic transport.

interaction between local moments due to flat bands. An approximate meeting ground of the two pathways lead to the Stoner criteria 43 , which argues that a high U and non-zero density of states (DOS) at Fermi energy is necessary for the unequal number of electrons with the two spins to be energetically favorable. Flat-band based mechanisms, or more generally the Stoner Criteria, for FM n −n , is thus supported only by the doped armchair edged graphene nano-ribbons(AGNR) on account of the flat or nearly flat bands which are located below and above the valence and conduction band edges and represent pi-bonds parallel to the ribbon edges. In zigzag edged graphene nano-ribbons(ZGNR), FM n −n has been proposed 44 to be possible with topological line defects owing to their non-bipartite nature.
Based on mean-field and ab-initio computation of spin resolved electronic structure, in this work we suggests an alternate approach to manipulate magnetism in graphene nano-ribbons(GNR), wherein, any GNR irrespective of its edge configuration, can be controllably as well as reversibly, turned into a FM semiconductor or metal and a half-metal, exclusively through spatially in-homogeneous biasing. Central to this approach is the robust emergence of FM n −n in the positively biased regions of the unit-cell, which we argue below to be a general property of in-homogeneously biased bipartite systems, arising primarily as a means to avoid increase of on-site Coulomb repulsion and kinetic energy while maximizing response to the external bias. This work points to a new class of mechanism for emergence of FM n −n order in bipartite systems in general at half filling, which is a clear departure from the body of work reported over the years on manipulation of magnetism. We demonstrate in the following the onset of FM n −n accompanied by lifting of spin degeneracy and consequent opening of window for ferromagnetic and half-metallic transport, in a representative variety of AGNRs, ZGNRs and also in a minimal model system.

Results and Discussion
We calculate spin polarized electronic structure of non-uniformly biased AGNR and ZGNR unit-cells shown in Fig. 1, first using the mean-field approximation of Hubbard model 31 as a function of bias potential V g and on-site Coulomb repulsion U, and then compare with density functional theory(DFT) 45 based first principles calculation.
To probe the nature of magnetic ordering we calculate the average n-n spin correlation as: where N s is the number of sites per unit-cell, nn i the number of n-n sites around the i-th site, and S i = 〈n i,σ 〉 − 〈n i,σ′ 〉, with 〈n i,σ(σ′) 〉 being the population of electron with spin-σ(σ′) at the i-th site due to the occupied states calculated using the mean-field approximation of Hubbard model. Positive and negative values of S plotted in Fig. 1 imply existence of FM and AFM or FeM ordering respectively. Existence of both thus imply spatial separation between FM and FeM ordering.
Negative spin correlation. For AGNRs, Fig. 1(a-f) imply rapid consolidation of AFM(FeM) ordering above U~2|t| with zero(positive) V g . For V g = 0 this is reminiscent of Mott transition 46 at half-filling (n = 1) in bipartite lattices. The trend that with increasing V g the transition from non-magnetic to the FeM ground state happens with increasing U, is similar to that observed in case of non-magnetic to AFM transition in bipartite lattices with increasing deviation from half-filling, and is understood in terms of the additional correlation required to dominate over the kinetic energy of the excess charges. The similarity in trend is expected since with non-zero V g the biased and unbiased regions both deviate locally from half-filling. With U > 0 at v g = 0, ZGNRs expectedly The degree of influence of width and location of biased region on the range of V g and U in which FM n −n would emerge, differs from ZGNRs to AGNRs. The difference can be attributed to the presence of weak but non-zero n-n FeM order in ZGNRs away from the edges. Wider biased region imply weaker localization, and thereby, weak spin separation. These results imply significance of localization and inherent magnetic order on emergence of FM n −n .
show n-n FeM ordering and AFM ordering globally between the two substructures. With increasing V g , quenching of magnetic ordering in ZGNRs below an increasing threshold of U can be understood as the dominance of positive bias over on-site Coulomb correlation, leading to occupation of biased sites by electrons with both the spins.
Positive spin correlation. Emergence of FM n −n ordering is marked by the positive spin correlation [ Fig. 1] and the associated lifting of spin degeneracy of the band-structures [ Fig. 2] over a range of V g with U moderate and higher. FM ordering quenches rapidly in AGNRs [ Fig. 1(a-f)] as the biased region moves away from the edges or are widened. Similarly in ZGNRs, positive spin correlation is much prominent if the biased region cover zigzag chains of carbon atoms parallel to the edges. Notably, an FM phase of generalized Nagaoka type is known to occur in the mean-field phase diagram of cubic bipartite lattice 31 at deviation from half-filling. Although Nagaoka may not be feasible 31 in three coordinated bipartite lattices, it is beyond the scope of this work to comment on whether it will be effective on the background of increased correlation due to confinement. However, the trend that the onset of the FM n −n order is more prominent if the biased region is narrow and located closer to the ribbon edges, clearly suggest that localization of electrons, and consequently the enhanced Coulomb correlation, are likely the key associated factors leading to FM n −n . Notably, in case of ZGNRs, if the biased sites cover zigzag(transverse) C-C bonds then the spin at those FM n −n ordered sites would prefer to be FM(AFM) ordered with sites at both the edges, mediated by the inherent n-n FeM order prevalent outside the biased region. Thus in case of biased zigzag sites, mixing of FM n −n ordered state with the dominant spin, and the nearest localized edge state, can stabilize both, leading to relative ease in occurrence of FM n −n order compared to that in case of biased transverse sites. This is consistent with less positive spin-correlation [ Fig. 1(j-l)] in case of biased C-C transverse bonds.
Band-structure. FM n −n order at the positively biased sites inherently implies lifting of spin-degeneracy since it consolidates one of the spins in the biased region. Emergence of the FM-semiconducting, FM-metallic and half-metallic phases are thus naturally expected as consequences of FM n −n . However, for AGNRs, onset of FM n −n order is preceded by shrinking of band-gap and direct to indirect transition. With increasing V g , the bands representing the bonding and anti-bonding orbitals of the biased region will come down in energy, while their counterparts associated with the edge farthest from the biased region will have least change in energies, resulting into a net reduction in band gap, as observed. Thus the conduction(valence) band edges can be spatially located at the positively(zero) biased regions of the unit-cell. In such a scenario a direct to indirect transition is expected, since the conduction and valence band edges, which are at different V g , tend to have same total energies as the band-gap shrinks, implying that their kinetic energies must be different. With increasing V g , further lowering of the conduction band edge leads to complete closure of the indirect band gap. The AGNRs are thus most likely to be normal or FM-metal upon emergence of FM n −n order [ Fig. 2(a,c)], although a small window for half-metallic transport can open up serendipitously. As the biased region moves towards the bulk(middle) of ribbons, lowering of kinetic energy increasingly compensates Coulomb correlation, leading largely to non-magnetic ground state and spin-degenerate band structure [ Fig. 2 In ZGNRs [ Fig. 2(d-g)] the effect of FM n −n ordering is particularly rich, since the FM ordering contests the inherent AFM ordering of the edge states by favoring same spin at the two edges, similar(opposite) to that of the FM n −n order if the biased region covers C-C zigzag(transverse) bonds. As shown schematically in Fig. 3, the contest supports the edge-states with spin-1 at edge-1, but suppresses the edge-state with spin-2 at edge-2. The consequent systematic closure of gap [ Fig. 2(f,g)] exclusively for spin-2, [Fig. 3] can be understood as a result of effective increase in on-site energy for electrons with spin-2 at edge-2 due to accumulation of electrons with spin-1 at that edge on account of FM n −n . Expedient to recall, equal probability of finding an electron at both the edges, implying same effective on-site term at both the edges, leads to closure of band-gap, as happens for both the spins in the absence of Coulomb correlation. Indeed, the reduction in n-n FeM order near edge-2 due to the competing magnetic orders, implies reduction of the effective Coulomb correlation near edge-2 where , the inter-edge AFM ordering clearly evolves into inter-edge FM ordering, leading to evolution of the valence band edge for spin-2 into a partially occupied dispersive band, which offers a robust window for half-metallic transport. This reiterating that bias coverage of zigzag C atoms is more effective than that of biased transverse C atoms in emergence of FM n −n order. However, the systematic reduction of band-gap for one of the spins, as argued above, happens in ZGNRs irrespective of whether the biased region covers C-C zigzag or transverse bonds, paving the way for robut bias controlled FM-semiconducting and half-metallic transport. Thus the same is offered by more realistic wider biased regions as well, at moderate U and low V g , which might be technologically relevant.
Notably, the emergence of the half-metallic window in this case is rooted at the emergence of the FM n −n order, and is thereby qualitatively different from previous proposals of bias induced opening of half-metallic window in ZGNRs where the sub-lattice asymmetry due to transverse electric field 13 or defects 24,25 , lead to FM ordering beyond nearest neighbours in a ferrimagnetic ground state. Furthermore, the robust emergence of the FM n −n order in the positively biased regions of the unit-cell, irrespective of location of the biased region with respect to the edges, and configuration of the ribbon edges, namely, zigzag or armchair or mix of both, clearly imply no explicit pre-requisition of flat or nearly flat bands at Fermi energy for the observed FM n −n order to emerge. Thus the spatially non-uniform bias induced FM n −n order proposed in this work is fundamentally different than the bias induced FM order 42 and half-metallic transport 13,26 , reported so far.
From first principles. Band-structure calculated from first principles using density functional theory (DFT) with sawtooth potential 45 applied in the transverse direction [ Fig. 4(a)] of ZGNR to resemble biased region akin to those considered in Fig. 1, has similar trend as [ Fig. 2(f)], marked by lifting of spin-degeneracy and reduction in band-gap for one of the spins [Fig. 4] within a range of ramp potential. This qualitative agreement between mean-field and DFT results is an important validation of the former, which only incorporates on-site Coulomb correlation. Furthermore, the reported agreement of DMRG, QMC and ED results 42,47 , with DFT in rationalizing FM ground state in doped AGNRs at moderate U, implies the confidence that our mean-field results will also be valid with improved consideration of correlation. Structural relaxation 48 using forces derived from DFT ground states appear to indicate that the ribbons should maintain their structure intact in the range of ramp potential for which the FM n −n order exists. However, in addition to the FM n −n ordering, we also find highly dispersive free electron like bands around Fermi energy [ Fig. 4(b)], whose origin is traced to accumulation of space charge between the periodic images of the ribbon. Upon emergence of FM n −n order in the ribbon, the space charge also acquires a net non zero magnetic moment.

Minimal model & mechanism.
To check the generality of our results beyond the three coordinated networks we resort to a minimal unit-cell, which can exhibit the FM n −n order if it exists. We choose a unit-cell of four consecutive sites, of which two neighboring sites are biased [inset, Fig. 5(b)]. Spin-correlation between the two biased sites [inset, Fig. 6] as a function of V g and U has similar generic features as in Fig. 1. thus hinting at the generality of the n-n FM order as a property of non-uniformly biased bipartite systems. Notably, if we do not consider non-zero crystal momentum, then the positive spin-correlation does not exist, although the trend of spin-correlation as a function of V g and U obtained using the mean-field Hubbard model, agrees qualitatively with that obtained using exact-diagonalization(ED). Itinerant electrons described by dispersive bands at Fermi energy are thus likely to play important role in manifestation of FM n −n order, as suggested by Fig. 2.   The evolution of spin polarized charge densities [ Fig. 5(a,c,e)] with increasing V g suggests spin-separation between biased and unbiased region similar to that observed in ribbons. Such a spin-separation can possibly be rationalized by noting that increased occupancy of the biased sites by electrons with both the spins, upon increased V g , would in turn increase potential energy due to on-site Coulomb repulsion. Thus to keep potential  (e-h). The light red and blue backgrounds used in (a) to (h) represent FeM and FM ordering respectively between the biased sites. The qualitative similarity of the energetics responsible for emergence of the FM n −n order, imply that the mean-field results are valid beyond the mean-field approximation. energy low, each of the biased sites will prefer to be dominated by electrons with one of the two spins. However if two such neighboring sites have opposite spins, as would be in a FeM ground state, then the wave-functions for both the spins will have rapid variation from site to site, leading to high kinetic energy. Instead, if wave-functions of one type of spin spans the biased sites more than those with the other spin, as evident from the charge densities [ Fig. 5(a,c,e)], then the wave-functions can be smoother than their FeM counterparts, implying lesser kinetic energy while allowing lower on-site Coulomb repulsion as well. The evolution of charge densities [ Fig. 5(a,c,e)] clearly implies electrons to be more itinerant for one spin than for the other upon emergence of FM n −n , leading to lifting of spin degeneracy [ Fig. 5(c,d)] akin to that in [ Fig. 2(a,c)]. Upon further increase of V g [ Fig. 6(b-d)] the lowering of potential energy due to V g dominates over the increase of potential energy due to on-site Coulomb repulsion, leading to occupancy of bias sites by both spins, resulting into return of non-magnetic ground state and spin degenerate band-structure [ Fig. 5(e,f)]. Thus with higher U a higher V g is required for the FM n −n order to quench, which is also consistent with the trend observed in Fig. 1.
To quantitatively justify the mechanism anticipated above, we partition the total energy [ Fig. 6(a)] of the unconstrained(UC) ground state into kinetic energy and potential energies due to on-site Coulomb repulsion and applied bias, and compare those with their counterparts in a non-magnetic(NM) and FeM ground states. NM ground state is obtained by assigning same charge density for both the spins. For the FeM ground state, charge density for one of the spins is assigned to be the mirror image of that for the other spin about the centre of the unit-cell. We choose the energetics of the NM ground state to be the reference. Figure 6(c,d) suggests that the UC ground state with FM n −n order has lower energy than NM and FeM ground states, owing initially to lowering of potential energy due to gate bias, but subsequently and primarily due to lowering of on-site Coulomb repulsion facilitated by separation of spin between positively biased and unbiased sites. Notably, although the FeM ground state has a lower kinetic energy than the UC ground state, the degree of localization at the positively biased site offered by the former is much lower than that due to the latter. Thus an FeM ordered state with the same degree of localization as offered by the UC ground state must have higher kinetic energy than the latter. Thus the observed FM n −n order is a result of attempts to minimize on-site Coulomb repulsion and kinetic energy in conjunction with each other, while maximizing localization of electrons at the positively biases sites.
In addition to agreement with DFT results, to demonstrate further that the observed FM n −n order is not limited by the mean-field approximation, we resort to simple analytical models for Wannier functions(WF), to represent the ground state of the four site unit-cell at half-filling. WFs are linear combination of wave-functions and can be chosen to be real and localized largely within a unit-cell in one dimension. Thus in place of 2N k wave-functions for each spin, N k being the number of allowed wave-vectors in the first Brillouin zone, we can choose two WFs for each spin to represent four electrons. We approximate WFs to be non-zero only within a unit-cell. Such approximate WFs describing the non-magnetic(NM) ground state can be of the following general form: which are constrained to be orthogonal to each other. Similarly, two orthogonal WFs to describe the FeM ground state can be approximated as: where the |φ i,↑ | 2 is mirror image of |φ i,↓ | 2 with respect to the middle of the unit-cell. Finally, orthogonal WFs with FM n −n order can be approximated as: The number of independent variables chosen to express the WFs are determined by the symmetry of the spin densities [ Fig. 5(a,c,e)] and orthogonality of the states. Total energies of ground states are calculated in the basis of the approximate WFs within the Hubbard model without mean-field approximation. For each class of WFs, ground state is obtained by finding the global minima of total energy using the cylindrical algebraic decomposition algorithm 49 . Kinetic energy and potential energies due to Coulomb repulsion and external bias are estimated using the WFs corresponding to the ground states. Notably, Fig. 6(e-h) implies emergence of FM n −n order in exactly the same pathway as suggested within the mean-field approximation of Hubbard model [ Fig. 6(a-d)]. These agreements indeed confirms the induced FM n −n order to be a genuine attribute of in-homogeneously ScIENtIfIc REPORTs | 7: 17094 | DOI:10.1038/s41598-017-17091-w biased bipartite systems and should be observed beyond mean-field computation, as also suggested by DFT. Thus the FM n −n order proposed in this work indeed appears to be a general property of non-uniformly biased bipartite systems. Although it remains to be specifically tested, available studies 44 suggest that similar emergence of FM n −n order may be even easier in non-bipartite systems.

Conclusion
To summarize, in this work we computationally propose the possibility of bias driven nearest neighbour(n-n) ferromagnetic(FM) order and half-metallic transport as a result of interplay of localization and Coulomb correlation. Results suggest emergence of n-n FM order in positively biased regions of non-uniformly biased bipartite systems, as a generic outcome of efforts to minimize Coulomb repulsion with minimal loss of the itinerant nature of electrons, while maximizing localization of electrons in the positively biased sites. The n-n FM order is computationally demonstrated here to exist in graphene nano-ribbons irrespective of their edge configurations, as well as in a minimal model system, which emphasizes the generality of the observed localization induced n-n FM order at half filling. The associated lifting of spin-degeneracy leads to exotic metallic (normal, ferromagnetic and half-metallic) phases. In armchair edged semiconducting ribbons the metallic phases are presided by direct to indirect transition, while in zigzag edged graphene ribbons, their inherent inter-edge AFM order contests the bias driven FM n −n order, leading to systematic closer of gap for one of the spins resulting into robust window for FM semiconducting and half-metallic transport. These results are expected to encourage a conceptually new pathway for voltage controlled opening of windows for half-metallic transport in two dimensional systems in general. In view of recent advancements [50][51][52] , in implementing gates at sub-micron length-scale, these results might encourage implementation of the proposed pathway for voltage controlled systematic opening of windows for half-metallic transport in ribbons and sheets of graphene and other two dimensional systems.

Method
To compute spin-polarized electronic structure in the realistic AGNR and ZGNR unit-cells, we resort to the mean-field approximation of Hubbard model within the n-n tight-binding framework: 〈n i,σ′ 〉 being the population of electron with spin-σ at the i-th site due to the occupied states. Eqn. (2) implies a self-consistent computation of electronic structure as a function of the strength of on-site Coulomb repulsion (U) 53,54 , between opposite spins, and the gate bias V { } i g , given the lowering of energy due to hopping between n-n sites represented by t = −2.7 eV [55][56][57] . Coulomb potential (V q ) due to electrons at nearest neighboring sites and beyond is calculated using the Ewald summation scheme 58 . We compare our mean-field Hubbard model based results with their counterparts obtained from first principles using density functional theory (DFT). We use a plane-wave based implementation 45 of DFT, wherein, we have used a gradient corrected functional 59 of density to approximately estimate the exchange-correlation contribution to total energy. Standard computational techniques have been used in this work to calculate electronic structure of systems detailed in the figures.