ZZ Polynomials for Isomers of (5,6)-Fullerenes C n with n = 20–50

: A compilation of ZZ polynomials (aka Zhang–Zhang polynomials or Clar covering polynomials) for all isomers of small (5,6)-fullerenes C n with n = 20–50 is presented. The ZZ polynomials concisely summarize the most important topological invariants of the fullerene isomers: the number of Kekulé structures K , the Clar number Cl , the ﬁrst Herndon number h 1 , the total number of Clar covers C , and the number of Clar structures. The presented results should be useful as benchmark data for designing algorithms and computer programs aiming at topological analysis of fullerenes and at generation of resonance structures for valence-bond quantum-chemical calculations.


Introduction
A (5,6)-fullerene is a polyhedral carbon cage with n carbon atoms arranged in 12 pentagonal and n 2 − 10 hexagonal faces [1]. The smallest (5,6)-fullerene consists of 20 carbon atoms and contains only pentagonal faces. The next smallest (5,6)-fullerenes are C 24 and C 26 , both possessing only a single isomer of point group symmetry D 6d and D 3h , respectively. All higher (5,6)-fullerenes C n with n ≤ 28 possess multiple isomers of various symmetry. The problem of generating all these isomers was solved first by introducing the ring spiral algorithm [2,3], valid for all the isomers of the (5,6)-fullerenes C n up to at least n = 200 [4], and later by the top-down approach of Brinkmann and Dress [5,6], valid in general cases. An invaluable compendium of useful information about all the isomers of (5,6)-fullerenes C n with n = 20-50 and so-called isolated pentagon rule (IPR) isomers of C n with n = 60-100 was compiled by Fowler and Manolopoulos in a form of a fullerene atlas [3], where the definition of each isomer is given as a sequence of 12 numbers denoting the positions of pentagons in the ring spiral. We follow the isomer labeling convention introduced by Fowler and Manolopoulos [3] also in the current work. Molecular structures in a XYZ format of all the isomers appearing in the Fowler and Manolopoulos atlas can be downloaded from "Fullerene Structure Library" [7]. (The reader should be warned that, for some mysterious reasons, geometries for the following pairs of isomers: (170,196) of C 48 and (44,178), (157,211), (27,59), and (115,170) of C 50 are switched in this library.) The remaining isomers can be conveniently generated using Fullerene, a program for the topological analysis of fullerenes written by Schwerdtfeger, Wirz, and Avery [8]. A compendium of various physical and chemical properties for (5,6)-fullerene isomers of C n with n = 20-42 were compiled previously by Małolepsza and collaborators [9,10], but this effort was discontinued once it was realized that for larger fullerene cages the differences between the studied properties are going to be minuscule, as they correlate mainly with local curvature of the fullerene cage rather than with the global topology of the carbon-carbon adjacency graph [11]. A wide-scope review summarizing many aspects of fullerene topology relevant in the context of the current report was given recently by Schwerdtfeger, Wirz, and Avery [12].
Various topological invariants of fullerene graphs received considerable attention in the literature . Most studies focused on the determination of the Clar number of fullerenes and on the computation of their number of Kekulé structures. Clearly, such a scope was motivated by practical considerations. Soon after the experimental discovery of C 60 , it was hypothesized [55] that the pronounced stability of the I h isomer of C 60 could be explained by a very large number of Kekulé structures that can be constructed for that isomer. This hypothesis was soon disproved, when Austin and collaborators discovered [1] that there exist 20 isomers of C 60 surpassing the count of 12,500 Kekulé structures for icosahedral C 60 (it might be relevant to mention here that only 158 of these Kekulé structures are symmetry distinct, i.e., not related by any I h point group symmetry operation [18]) and demonstrated the lack of obvious correlation between the raw Kekulé number K and other, quantum-chemical descriptors of its pronounced stability. We mention in passing that the isomer of C 60 with maximal K has 16,501 Kekulé structures, which is considerably larger (by some 30%) than for the I h isomer. Interestingly, a recent accurate study of the thermodynamic stability for the isomers of C 60 shows that the vast majority of the isomers with K > 12,500 belong to the most thermodynamically unstable isomers of C 60 [48]. Only relatively recently was it understood that the correct perspective comes not only from looking at the number of Kekulé structures, but also at their Clar numbers, i.e., the maximal number of aromatic Clar sextets [56] that can be accommodated by the fullerene graph. Zhang and collaborators demonstrated [57] that the icosahedral C 60 indeed has the highest Kekulé count among the isomers of C 60 with the largest Clar number, Cl = 8. There exists 18 of such isomers and the second highest Kekulé count among them is 11,259, about 10% lower than for the I h isomer. One should not, however, overemphasize this results, as most likely the pronounced stability of the icosahedral C 60 has not only thermodynamic but also kinetic provenance [58].
It would be also interesting to correlate topological indices of fullerene isomers with their thermodynamic stabilities for fullerenes other than C 60 . The main problem for such a manifesto is the lack of data allowing for such comparisons. The main motivation for the current study is filling this gap by compiling a collection of topological indices for all the isomers of small (5,6)-fullerenes with n ≤ 50. The current work can be considered as an extention of the tabulation of matching polynomials given by Balasubramanian [14] We would be happy to extend this compilation also to larger fullerenes (particularly to the isomers of C 60 , which occupy a pronounced position in practical considerations), but the spatial extent of such a tabulation would exceed any sensible length advisable for a scientific paper. The topological indices are given in the form of ZZ polynomials. This choice is rather clear, as ZZ polynomials are probably the most concise and robust form of presenting such invariants available in the literature, containing information about the number of Kekulé structures K, the Clar number Cl, the first Herndon number h 1 , the total number of Clar covers C, and the number of Clar structures for each of the isomers, in addition to the number of Clar covers of each order. Detailed information about ZZ polynomials, their structure, and their way of determination are given in the next section.

ZZ Polynomials
ZZ polynomials (aka Zhang-Zhang polynomials or Clar covering polynomials) were introduced to the field of chemical graph theory about 25 years ago by two Chinese mathematicians, Fuji Zhang, and Heping Zhang [59][60][61][62]. Formally speaking, a ZZ polynomial ZZ(B, x) of some benzenoid structure B is a generating function for the sequence of the numbers of Clar covers of B of each order. Since this definition is not widely known among chemists, let us introduce here the concept of a ZZ polynomial on a simple example of benzo[e]pyrene. Figure 1 shows (in black) the molecular structure of benzo[e]pyrene together with the entire collection of Clar covers (in gray) that can be constructed for this molecule. A Clar cover is a generalized resonance structure, in which the tetravalent character of each carbon atom has been satisfied by distributing a certain number (say m) of double bonds and a certain number (say k) of aromatic Clar sextets. Since each double bond involves two carbon atoms and each aromatic sextets involves six carbon atoms, and since each carbon atom can be involved in only one of those structures, we have a natural connection between the number n of carbon atoms in B and the non-negative integers m and k, given by 2m + 6k = n. This condition shows that 0 ≤ k ≤ n 6 and that m = n 2 − 3k. It often happens in benzenoid hydrocarbons that the natural upper bound for k given by n 6 is not achieved, as it is not possible to arrange the double bonds in a way that is compatible with n 6 aromatic sextets. Benzo[e]pyrene consists of n = 20 carbon atoms and the number of aromatic sextets that can be accommodated in this molecule is 0 ≤ k ≤ 3. All of these possibilities are realized. Figure 1 shows that there exists exactly one Clar cover (usually referred to as the Clar structure or Clar formula) of benzo[e]pyrene with k = 3 aromatic sextets (depicted in a blue frame in Figure 1). Therefore, the Clar number Cl of this molecule is 3. Similarly, we can construct seven distinct Clar covers of order 2 (in the green frame in Figure 1), each of them comprising two aromatic sextets and four double bonds. The corresponding numbers of distinct Clar covers of order 1 and 0 are 16 and 11, respectively. Clearly, the Clar covers of order 0 (in the orange frame in Figure 1 Clearly, the ZZ polynomial of benzo[e]pyrene decodes its most important topological invariants: the number of Kekulé structures K = c 0 , the Clar number Cl = degree(ZZ(benzo[e]pyrene, x)), the first Herndon number h 1 = c 1 [63], the total number of Clar covers C = ZZ(benzo[e]pyrene, 1), and the number of Clar structures is given by c Cl . The most attractive feature of ZZ polynomials is the robustness with which they can be determined. Computing a single coefficient in the ZZ polynomial or determination of the Clar number of a given benzenoid is a computationally complex problem. However, determination of the entire ZZ polynomial is much simpler owing to the convenient recursive properties it obeys. Zhang and Zhang in the original paper [59] derived a number of decompositions (see Theorems 3-6 of [59]) allowing for computing the ZZ polynomial of a given benzenoid B as a weighted sum of its substructures. Building on this principle, it is possible to design a recursive algorithm, which performs such a decomposition multiple number of times and computes the final ZZ polynomials from the ZZ polynomials of the nodes of the recursive decomposition tree. Such an algorithm was first proposed by Gutman and collaborators [64] and the details and an actual robust implementation were provided by our group [65,66]. The original program, written in Fortran90, was soon accompanied by a graphical engine ZZDecomposer [67][68][69] allowing for visualizing in real time the recursive decomposition pathways and allowing for discovering closed-form ZZ polynomial formulas for the whole families of isostructural benzenoids [70]. At the moment, such closed-form ZZ polynomial formulas are known for almost all of the families of basic benzenoids [64,[71][72][73][74][75][76][77][78][79][80][81][82] with the exception of hexagonal graphene flakes and oblate rectangles, which are the objects of current intensive research activity [83][84][85]. Another interesting property of ZZ polynomials is their equivalence to cube polynomials [42,86,87] and certain tiling polynomials [88].
The concept of a ZZ polynomial almost immediately can be generalized from benzenoid hydrocarbons to fullerenes. The only limitation to be imposed on the algorithm constructing the Clar covers of fullerenes is that the aromatic sextet-obviously-cannot be placed in any of the fullerene's pentagons; the remaining rules for determination of ZZ polynomial coefficients are exactly the same as for benzenoid hydrocarbons. The results reported in the next section are computed using a stand-alone Fortran90 code [67] with executables included in every distribution of ZZDecomposer [68,69]. The source code of the program can be obtained from the authors upon request. Few of the computed ZZ polynomials have been verified by pencil-and-paper calculations to make sure that no programming errors are present in the used subroutines.

List of ZZ Polynomials for Fullerene Isomers
The computed ZZ polynomials for all the isomers of (5,6)-fullerenes C 20 -C 50 are presented in Table 1. For each of the isomers, in addition to its point group symmetry designation, we give its two distinct definitions. The first definition is given in the form of an isomer number identical to the designation presented in "An atlas of fullerenes" compiled by Fowler and Manolopoulos [3]. This reference also gives the ring spiral pentagon sequence for each isomer, which can be used to generate its geometry in the XYZ format with the Fullerene program [8]. The isomer number can be also used to extract the XYZ geometry of the isomer from the "Fullerene Structure Library" [7]. The second definition is given in the form of a Schlegel diagram, obtained by taking the optimized geometry of each isomer, projecting it in the surface of a unit sphere and subsequently projecting the thus produced points of the sphere to a 2D plane by a stereographic projection. The geometrical structures of all the fullerene isomers have been optimized using density-functional tight-binding code, like in our previous publications [9][10][11]89] on this topic.
It is interesting to note that all the ZZ polynomials presented here are distinct from each other. This is in clear contrast to single topological indices (e.g., the Clar number Cl or Kekulé count K) of these isomers, which often have the same value for different isomers.
Let us explain on the example of the isomer 22 of C 50 -denoted further for convenience as C 50 : 22-how various topological invariants of the isomers can be extracted from its ZZ polynomial. According to Table 1, the ZZ polynomial of C 50 : 22 is given by the following expression: ZZ(C 50 : 22, x) = 2541 + 4820x + 3478x 2 + 1222x 3 + 216x 4 + 16x 5 (2) The following topological invariants can be readily extracted from this expression: • Kekulé count K is equal to the coefficient of x 0 , so here we have K = 2541. Note that one can alternatively evaluate the ZZ polynomial at x = 0 to obtain the same value. • Clar number Cl is equal to the degree of the ZZ polynomial, so here we have Cl = 5.

•
The total number C of Clar covers is equal to the sum of all the coefficients in the ZZ polynomial. C is most conveniently computed by evaluating the ZZ polynomial at x = 1. For C 50 : 22, we have C = ZZ(C 50 : 22, 1) = 12293.

•
The number of Clar formulas, i.e., the number of Clar covers with the maximal number Cl of aromatic sextets, is equal to the coefficient of x Cl , which for C 50 : 22 is equal to 16.

•
The first Herndon number is equal to the coefficient of x 1 , which for C 50 : 22 is equal to 4820. Table 1. Compilation of ZZ polynomials for all the isomers of small (5,6)-fullerenes C n with n = 20-50. The columns specify: the fullerene type, the isomer number (following the convention introduced in [3]), point group symmetry, Schlegel diagram, The Kekulé count K, the ZZ polynomial, and the total number C of Clar covers of a given isomer.

Discussion
Our main motivation here is to provide the chemical graph theory community with a compilation of reference data for graph-theoretical invariants of all isomers of small fullerenes C 20 -C 50 . However, having access to such a collection of data, we cannot resist-stimulated by a request of an anonymous referee-to show the correlations between the computed invariants and energetic stability of each isomer in spirit of numerous publications devoted to this topic [1,13,16,52,59,[90][91][92][93][94][95][96][97]. To this end, we have optimized the geometry of each fullerene isomer using a density-functional tight-binding (DFTB) method [98] following the methodology of our previous papers [9][10][11]. To facilitate the comparisons between different fullerenes, each of the optimized energies is divided by the number of carbon atoms in this structure, obtaining a descriptor of thermodynamic stability of each isomer in a form of energy per carbon atom. The structure with the lowest energy corresponds to the most stable isomer of a given fullerene. The comparison between the energies (per atom) of each isomer and their Kekulé count K and Clar count C are given in Figures 2 and 3, respectively. Briefly summarizing these results, one can say that in general K and C are very loosely correlated with the computed energies. The isomers with the highest values of K for each fullerene usually correspond to thermodynamically unstable structures characterized often by the highest energies. Organic chemists often relate the stability of a given species to the number of resonance structures that can be drawn for it; structures with more resonance forms are thought to be more stable than those with less resonance forms. Apparently, this logic cannot be extended to the fullerene isomers, in agreement with previous observations made for C 60 [1,48], as mentioned in Section 1. (Our results can be compared also with previous study of Manoharan and collaborators [99].) The most stable structure of a given isomer corresponds usually to an intermediate value of K. Similar observations can be made for the correlation between the computed energies of each isomer and their Clar count C. The main difference between the plots shown in Figures 2 and 3 concerns the vertical distribution of the circles; for K, they are negatively correlated with energies (i.e., a higher value of K corresponds to a higher, less-stable energy), while, for C, there seem to be no particular correlation of such a type.
We have mentioned previously in Section 1 that Zhang and collaborators demonstrated [57] that the most stable, icosahedral isomer of C 60 maximizes the Kekulé count among the isomers of C 60 with the largest Clar number, Cl = 8. It is indeed an interesting observation worth testing also for other fullerenes. Here, in Figure 4, we have performed an analogous analysis for two smaller fullerenes, C 36 and C 50 . Their DFTB energies per atom (multiplied by −1 and expressed in atomic units) are plotted as vertical bars with the position of each bar determined by the Clar number Cl and Kekulé count of a given isomer. To facilitate the comparisons, the three most stable isomers of each fullerene are designated with numbers 1 , 2 , and 3 . The most important observation concerning the presented data are that the performed analysis for C 36 and C 50 does not support the conclusion drawn by Zhang and collaborators [57] for C 60 . The two most stable isomers of C 36 have the smallest Clar number (equal to 2) among all the isomers of C 36 50 , the most stable isomer indeed has the highest Clar number (equal to 5) among the isomers of C 50 , but its Kekulé count has rather an intermediate value. The next two most stable isomers of C 50 , with Clar numbers of 4 and 5, respectively, are characterized even by a smaller Kekulé count. The isomer maximizing the Kekulé count among the isomers of C 50 with the largest Clar number, Cl = 5, actually belongs to the most unstable isomers of C 50 . All these data suggest that the correlation between stability and the maximal Kekulé count among isomers with the highest Clar number observed for C 60 by Zhang and collaborators has been rather accidental and cannot be immediately generalized to other fullerenes, but a definitive conclusion of that kind would require repeating the performed here analysis with more accurate quantum chemical technique than the inherently approximate DFTB method used here.

Conclusions
We have presented a short introduction to the theory of ZZ polynomials (aka Zhang-Zhang polynomials or Clar covering polynomials) accompanied by a compilation of ZZ polynomials for all the isomers of small (5,6)-fullerenes C 20 -C 50 . The ZZ polynomials concisely summarize the most important topological invariants of the fullerene isomers: the number of Kekulé structures K, the Clar number Cl, the first Herndon number h 1 , the total number of Clar covers C, and the number of Clar structures. The results are presented in a form of an extended table, where every isomer is identified by giving a reference number corresponding to the fullerene atlas [3] and by its Schlegel diagram. We believe that the presented results should be useful as benchmark data for designing algorithms and computer programs aiming at topological analysis of fullerenes and upon the generation of resonance structures for valence-bond quantum-chemical calculations.