Electronic Properties of Hexagonal Graphene Quantum Rings from TAO-DFT

The reliable prediction of electronic properties associated with graphene nanosystems can be challenging for conventional electronic structure methods, such as Kohn–Sham (KS) density functional theory (DFT), due to the presence of strong static correlation effects in these systems. To address this challenge, TAO (thermally assisted occupation) DFT has been recently proposed. In the present study, we employ TAO-DFT to predict the electronic properties of n-HGQRs (i.e., the hexagonal graphene quantum rings consisting of n aromatic rings fused together at each side). From TAO-DFT, the ground states of n-HGQRs are singlets for all the cases investigated (n = 3–15). As the system size increases, there should be a transition from the nonradical to polyradical nature of ground-state n-HGQR. The latter should be intimately related to the localization of active TAO-orbitals at the inner and outer edges of n-HGQR, which increases with increasing system size.

Due to the improvement in graphene fabrication, GQRs have attracted increasing attention in recent years. Note that GQRs are closely related to GQDs, and can be regarded as GQDs with additional inner edges. It is expected that the inter-edge coupling of states and enhanced finite-size effects in GQRs could yield novel properties. Indeed, some interesting properties of GQRs, such as the Aharonov-Bohm (AB) effects [42], have been theoretically predicted [39] and experimentally observed [34]. Further theoretical studies have focused on the AB effects on the persistent currents in GQRs [44] as well as the interplay of the AB effects and Klein tunneling in GQRs [45]. In addition, other theoretical studies have revealed the importance of edge and geometry effects on the energy spectra of GQRs in a magnetic field [36,43].
Among various GQRs, in this study, we focus on the hexagonal GQRs, consisting of n aromatic rings fused together at each side, which are denoted as n-HGQRs. As shown in Figure 1, n-HGQR (with the chemical formula C 24(n−1) H 12(n−1) ), which contains N = 156(n − 1) electrons, has six-fold rotational symmetry. To date, there have been scarce experimental studies on n-HGQR, which may be due to the fact that it remains very difficult to realize the precise n-HGQRs with smooth edges. On the theoretical side, there have been only few studies on the properties of energy spectra associated with n-HGQRs in the presence of a magnetic field [46,47]. Owing to the presence of strong static correlation in the larger n-HGQRs (which will be shown and discussed later), conventional electronic structure methods may not be appropriate. For example, Kohn-Sham (KS) density functional theory (DFT) [48,49] using the common exchange-correlation (XC) energy functionals (e.g., the LDA (local density approximation), GGA (generalized-gradient approximation), hybrid as well as doublehybrid XC functionals) have been commonly applied to predict the properties of electronic systems without strong static correlation (i.e., electronic systems with a single-reference (SR) character or SR systems), while KS-DFT with the common XC functionals can provide unreliable predictions of the properties of electronic systems with strong static correlation (i.e., electronic systems with a multi-reference (MR) character or MR systems). On the other hand, although high-level ab-initio MR electronic structure methods [50][51][52] are expected to take care of strong static correlation, accurate MR methods are known to be computationally expensive, and, thus, inapplicable to nanosystems, such as the larger n-HGQRs.
Recently, thermally assisted occupation (TAO) DFT, which is a DFT using fractional TAO-orbital occupation numbers (obeying the Fermi-Dirac (FD) distribution with a fictitious temperature θ) has been developed to describe the properties of nanosystems with strong static correlation [53]. Unlike finite-temperature DFT (FT-DFT) [54] (i.e., a DFT for studying the thermal equilibrium properties associated with systems at finite electronic temperatures; FT-DFT has been powered by the Mermin-Kohn-Sham (MKS) equations (see Section III of Ref. [49]), and the resulting MKS theory [49,54] is reduced to KS-DFT at zero electronic temperature (absolute zero)), TAO-DFT is a DFT for studying the groundstate properties associated with systems at absolute zero (even with a nonzero fictitious temperature θ, i.e., a parameter closely dependent on the strength of the static correlation of an electronic system at absolute zero). In the framework of TAO-DFT, the static correlation energy associated with the ground state of an electronic system is treated explicitly (though approximately) by the entropy contribution, and the corresponding natural orbital occupation numbers (NOONs) are approximately given by the TAO-orbital occupation numbers (TOONs) [53]. For the limit case of θ = 0, TAO-DFT is reduced to KS-DFT. It is worth mentioning that TAO-DFT maintains a similar computational complexity as KS-DFT, while TAO-DFT can outperform KS-DFT for a diverse range of MR systems (even with the simplest LDA functional).
As a result, in this work, we use TAO-DFT to explore the electronic properties (e.g., the singlet-triplet gaps, fundamental gaps and vertical electron affinities/ionization potentials, symmetrized von Neumann entropy as well as the occupation numbers and visualization of active TAO-orbitals) of n-HGQRs (n = 3-15). One of the objectives of this study is to address whether there is a transition from the nonradical to polyradical nature of ground-state n-HGQR with increasing system size, which has never been reported, although such a transition has been previously found in hexagonal GQDs [68] (i.e., closely related electronic systems). The significance of this work is to address whether unconventional electronic structure methods, such as TAO-DFT, are necessary for a reliable prediction of electronic properties of n-HGQRs.

Computational Details
We carried out all numerical calculations with the 6-31G(d) basis set on the software of Q-Chem 4.4 [78]. Geometry optimizations (see Supplementary Information (SI)) as well as single-point energy calculations were carried out using TAO-LDA, which is TAO-DFT [53] with the LDA XC and θ-dependent functionals, wherein θ = 7 mhartree (i.e., an optimal system-independent θ value for TAO-LDA calculations on both SR and MR systems) [53] was adopted. For comparison purposes, we also present the results of KS-LDA (i.e., KS-DFT using the LDA XC functional or TAO-LDA with θ = 0).

Singlet-Triplet Gap
To obtain the ground state of n-HGQR, we resorted to the singlet-triplet (ST) gap. First, we carried out KS-LDA and TAO-LDA calculations for the lowest spin-unrestricted triplet energy (E UT ) and lowest spin-unrestricted singlet energy (E US ) of n-HGQR at the respective optimized geometries. Subsequently, we computed the ST gap (E ST ) of n-HGQR using As shown in Figure 2, TAO-LDA predicts that n-HGQRs (n = 3-15) (i.e., all the cases studied) should possess singlet ground states (see also Table S1 and Figure S1  To explore the possible reasons for discrepancies, we reveal the expectation values of the total spin-squared operator Ŝ 2 associated with the lowest singlet and lowest triplet states of n-HGQR, computed using spin-unrestricted KS-LDA (see Table 1). If the KS-LDA wavefunction has spin contamination (which is the artificial mixing of different spin-states) [68,71,79], the corresponding Ŝ 2 value will be larger than the exact Ŝ 2 value (i.e., 0 for the lowest singlet state and 2 for the lowest triplet state; e.q., refer to Section I in the SI of Ref. [79]). Consequently, the difference between the Ŝ 2 value of KS-LDA and the exact Ŝ 2 value is commonly used to assess the degree of spin contamination in the KS-LDA wavefunction. As shown, for the smaller n-HGQRs, the Ŝ 2 values of KS-LDA are rather close to the exact Ŝ 2 values, indicating that the lowest singlet and lowest triplet states of n-HGQR, obtained with spin-unrestricted KS-LDA, essentially have no spin contamination. Nevertheless, for the larger n-HGQRs (e.g., n ≥ 11), the Ŝ 2 values of KS-LDA are much larger than the exact Ŝ 2 values, indicating that the lowest singlet and lowest triplet states of n-HGQR, obtained with spin-unrestricted KS-LDA, are highly spin-contaminated. Therefore, the larger n-HGQRs (e.g., n ≥ 11) are expected to possess a pronounced MR character in the lowest singlet and lowest triplet states, and the unexpected increase in the ST gap of KS-LDA beyond 11-HGQR can be an artifact intimately correlated with spin contamination.
In addition, we also examine the consequences of spin contamination in spin-unrestricted KS-LDA/TAO-LDA calculations for the lowest singlet states of n-HGQRs. Owing to the constraint of spin symmetry in the lowest singlet state of n-HGQR, for the exact theory, the lowest singlet energy of n-HGQR computed using the spin-restricted formalism must be the same as that computed using the spin-unrestricted formalism [53,55,56]. Accordingly, we carried out KS-LDA and TAO-LDA calculations for the lowest spin-restricted singlet energy (E RS ) and lowest spin-unrestricted singlet energy (E US ) of n-HGQR at the corresponding optimized geometries. Subsequently, we computed the E UR value of n-HGQR using (2) As presented in Table 2, the E UR values of the larger n-HGQRs (e.g., n ≥ 11), calculated by KS-LDA, are very large (e.g., more than 1.6 kcal/mol), leading to unphysical spin-symmetry breaking effects in the respective spin-unrestricted KS-LDA calculations. By contrast, the E UR values of n-HGQRs (n = 3-15) (i.e., all the cases examined), calculated by TAO-LDA, remain essentially zero, yielding essentially no unphysical spin-symmetry breaking effects in the respective spin-unrestricted TAO-LDA calculations. In short, the ST gaps associated with the larger n-HGQRs (e.g., n ≥ 11), calculated by KS-LDA, can be seriously affected by spin contamination, greatly degrading the accuracy of KS-LDA in predicting the ST gaps and possibly other electronic properties of the larger n-HGQRs (e.g., n ≥ 11). As will be shown later, the larger n-HGQRs possess an increasing polyradical nature in the ground states (i.e., the lowest singlet states), wherein KS-LDA can provide unreliable electronic properties. Therefore, we only report the results of TAO-LDA hereafter.

Fundamental Gap and Vertical Electron Affinity/Ionization Potential
Another important electronic property is the fundamental gap, which is known to vanish for metals, but is positive for semiconductors and insulators [80]. Here, we employ spin-unrestricted TAO-LDA to calculate the fundamental gap (E g ) of ground-state n-HGQR (containing N electrons): which is the difference between the IP v (vertical ionization potential) and EA v (vertical electron affinity) of ground-state n-HGQR. To obtain the IP v /EA v of ground-state n-HGQR, we perform spin-unrestricted TAO-LDA calculations with one electron removed from/added to neutral n-HGQR at the ground-state geometry, and compute the following energy difference [55,56]: and where E N refers to the total energy associated with the N-electron system. As shown in Figure 3, with the increase in system size n, the E g of ground-state n-HGQR is monotonically decreasing, which is a consequence of the monotonically decreasing IP v and monotonically increasing EA v of ground-state n-HGQR (see, also, Table S2 in SI). Moreover, the fundamental gaps of ground-state n-HGQRs (n = 5-10) are in the range of 1-3 eV (i.e., for the ideal energy gaps), which shows great promise for the applications of ground-state n-HGQRs (n = 5-10) in nanophotonics.

Symmetrized von Neumann Entropy
To assess the strength of MR character associated with ground-state n-HGQR, we adopted spin-unrestricted TAO-LDA to calculate the S vN (symmetrized von Neumann entropy) associated with ground-state n-HGQR [55,56]: Here, f i,σ (a value in the range of 0-1) is the occupation number associated with the ith TAO-orbital of σ-spin (i.e., up-or down-spin), which is intimately correlated with the respective natural spin-orbital occupation number [53,55,56]. For an SR system, the f i,σ values associated with all TAO-spin-orbitals must be very close to either 0 or 1, yielding a very small value of S vN . However, for an electronic system with pronounced MR character, the f i,σ values associated with the active TAO-spin-orbitals (i.e., the TAO-spin-orbitals with significant fractional occupations (e.g., 0.1-0.9)) can strongly deviate from the values of 0 and 1; hence, the corresponding value of S vN can noticeably increase when the f i,σ values associated with the active TAO-spin-orbitals are closer to 1 2 , and/or the number of active TAO-spin-orbitals increases. As a consequence, the value of S vN associated with ground-state n-HGQR can be adopted to assess the strength of MR character associated with ground-state n-HGQR. It is worth mentioning that the strength of MR character associated with ground-state n-HGQR may also be explored by other measures [81,82], due to the availability of TOONs (i.e., approximate NOONs) in TAO-DFT.
As shown in Figure 4, the S vN of ground-state n-HGQR, which remains very small for a small value of n, is monotonically increasing with the increase in n-HGQR size (see Table S2 in SI as well). This implies that the smaller ground-state n-HGQR should exhibit an SR character (i.e., nonradical nature), and the larger ground-state n-HGQR should possess an MR nature whose strength generally increases with increasing system size n.  Finally, we briefly comment on the fractional occupation number weighted density (FOD) method [83], which has been recently proposed for the real-space measure of strong static correlation effects at absolute zero. In the FOD method, the special density ρ FOD is actually computed using an approximate TAO-DFT method [53] (i.e., TAO-DFT without the θ-dependent functional) at some fictitious temperature θ, instead of using an approximate MKS method [49,54] (i.e., the MKS theory with the XC free energy functional being approximated with the XC functional) at finite electronic temperature. This can be clearly seen from the distinctly different physical meanings of TAO-DFT [53] and the MKS theory [49,54] (which is reduced to KS-DFT at absolute zero). In other words, the fractional orbital occupation numbers adopted in the FOD method are approximate TOONs, and, hence, can be viewed as approximate NOONs. Accordingly, the N FOD value (obtained from the integration of ρ FOD over all space), which is similar to the aforementioned S vN value, can also be used to measure the strength of the MR character associated with the ground state of an electronic system at absolute zero (even with a very large fictitious temperature θ), rather than that associated with an electronic system in thermal equilibrium at a finite electronic temperature.

Active TAO-Orbital Occupation Numbers
As mentioned previously, the NOONs associated with the ground state of an electronic system are approximately given by the TOONs in TAO-DFT [53,55,56]. To provide further insight into the increase in S vN with increasing n-HGQR size and to understand the radical nature of the larger ground-state n-HGQR, we disclose the active TOONs of ground-state n-HGQR, obtained from spin-restricted TAO-LDA. For ground-state n-HGQR (containing N electrons), we denote the HOMO (highest occupied molecular orbital)/LUMO (lowest unoccupied molecular orbital) as the (N/2)th/(N/2 + 1)th TAO-orbital, and so forth. For brevity, we denote HOMO/LUMO as H/L. In addition, the TAO-orbitals with an orbital occupation number in the range of 0.2-1.8 (i.e., the TAO-spin-orbitals with an orbital occupation number in the range of 0.1-0.9) are viewed as the active TAO-orbitals.
As shown in Figure 5, the active TOONs of ground-state n-HGQR display peculiar six-member patterns (i.e., the active TOONs are arranged to form groups of six members each), which can be related to the six-fold rotational symmetry of n-HGQR. When the n-HGQR size is small (e.g., n < 7), all the TOONs are close to either 0 or 2. Therefore, the smaller ground-state n-HGQRs (e.g., n < 7) should possess a nonradical nature, showing consistency with the investigation of the other properties, such as the larger E ST and E g values and the smaller S vN values, associated with the relatively stable n-HGQRs. Nevertheless, as the n-HGQR size increases, the active TOONs get closer to 1, and/or the number of active TAO-orbitals increases. Accordingly, the larger ground-state n-HGQRs (e.g., n ≥ 7) should possess an increasingly polyradical nature, also showing consistency with the investigation of the other properties, such as the smaller E ST and E g values and the larger S vN values, associated with the relatively unstable n-HGQRs. As the system size grows, there is a transition from the nonradical to polyradical nature of ground-state n-HGQR, leading to the aforementioned increase in S vN .

Real-Space TAO-Orbital Representation
To further illustrate the transition from the nonradical to polyradical nature of groundstate n-HGQR in real space, we plot the visualization of active TAO-orbitals (H−5, H − 4, . . . , H − 1, H, L, L + 1, . . . , L + 4, and L + 5) of ground-state n-HGQRs (n = 3, 5, 7, and 9), calculated by spin-restricted TAO-LDA. As shown in Figures 6-9, when the n-HGQR size is small (e.g., n < 7), the active TAO-orbitals are delocalized over the entire n-HGQR. Nonetheless, as the n-HGQR size increases, the active TAO-orbitals have an increasing tendency to localize at the inner and outer edges of n-HGQR. Consequently, the active TAO-orbitals of the larger ground-state n-HGQR (e.g., n ≥ 7) are mainly localized at the two edges of n-HGQR.

Conclusions
In conclusion, the MR character of graphene nanosystems has posed a great challenge to conventional electronic structure methods. The challenge has been greatly reduced by TAO-DFT. In this study, we explored the electronic properties of n-HGQRs (n = 3-15) obtained with TAO-LDA. Since the larger n-HGQRs were found to have pronounced MR character in the lowest singlet and lowest triplet states, KS-DFT using the common XC functionals can be unreliable for predicting the electronic properties of n-HGQRs. For example, the lowest singlet and triplet states of the larger n-HGQRs (e.g., n ≥ 11), obtained with spin-unrestricted KS-LDA, were found to be highly spin-contaminated, leading to unphysical spin-symmetry breaking effects. Although accurate MR electronic structure methods are able to describe strong static correlation effects, these methods can be computationally intractable for the larger n-HGQRs. Consequently, it is reasonably justified that we employ TAO-LDA to study the electronic properties of n-HGQRs due to a good compromise between efficiency and accuracy.
From the TAO-LDA results, n-HGQRs (n = 3-15) are all ground-state singlets. With increasing n-HGQR size, the E ST , IP v , and E g values monotonically decrease, while the EA v and S vN values monotonically increase. With increasing system size, there is a transition from the nonradical to polyradical nature of ground-state n-HGQR. The increasing polyradical character associated with the larger ground-state n-HGQR (e.g., n ≥ 7) should be intimately related to the localization of active TAO-orbitals at the inner and outer edges of n-HGQR, which increases with increasing system size.
As for possible future work, we plan to explore how the radical nature emerges from the evolution between n-HGQRs and hexagonal GQDs (i.e., the coronene series) [68] to see the trend forming hexagonal GQDs.

Data Availability Statement:
The data that support the findings of this study are available from the authors upon reasonable request.