Spectral clustering of combinatorial fullerene isomers based on their facet graph structure

After Curl, Kroto and Smalley were awarded 1996 the Nobel Prize in chemistry, fullerenes have been subject of much research. One part of that research is the prediction of a fullerene's stability using topological descriptors. It was mainly done by considering the distribution of the twelve pentagonal facets on its surface, calculations mostly were performed on all isomers of $C_{40}, C_{60}$ and $C_{80}$. This paper suggests a novel method for the classification of combinatorial fullerene isomers using spectral graph theory. The classification presupposes an invariant scheme for the facets based on the Schlegel diagram. The main idea is to find clusters of isomers by analyzing their graph structure of hexagonal facets only. We also show that our classification scheme can serve as a formal stability criterion, which became evident from a comparison of our results with recent quantum chemical calculations. We apply our method to classify all isomers of $C_{60}$ and give an example of two different cospectral isomers of $C_{44}$. Calculations are done with MATLAB. The only input for our algorithm is the vector of positions of pentagons in the facet spiral. These vectors and Schlegel diagrams are generated with the software package Fullerene.


Introduction
In this paper, we consider a fullerene C n as a convex polytope modeling a closed three-dimensional carboncage with n atoms, cf. [1]. Each vertex is connected to exactly three other vertices, such that the facets are pentagons and hexagons only. Using the classical Euler relation and Eberhard's theorem [19, §13.3] for simple three-dimensional convex polytopes one can conclude that the number of pentagonal facets is always equal to 12 and n is even. In this case, the number of hexagons is m 6 = n 2 − 10, cf. [17]. Theoretically, C n exists for n = 20 and all even n ≥ 24, see [1]. In the sequel we call such a n feasible. However, up to now Figure 1: IPR-Isomers of C 60 (Buckminster fullerene) and C 70 . Courtesy of Max von Delius only a few of them (some isomers of C n with n = 60, 70, 76, 78, 80, 82, 84) are separated (i.e., are chemically stable and can be synthesized in considerable mass quantities), cf. [21,27,30,38]. Another difficult problem is the combinatorial constructive enumeration of all isomers of C n , see [8,12,35] and references in [16,Section 3]. For instance, a set of four operations (including that of Endo-Kroto) enables the construction of all fullerenes for each even n starting with n = 24 (barrel), see [13,16]. Other examples for such operations are the generalized Stone-Wales operation (cf. [2], which we investigate in a forthcoming paper [5] in more detail) and the buckygen (introduced 2012 in [9]) which is up to now the fastest algorithm to generate all C n -isomers.
Combinatorial isomer is a class of combinatorially equivalent polytopes. Two polytopes P and P are called combinatorially equivalent if there exists a one-to-one mapping between the lattice of all facets of P and P that is inclusion-preserving [19, p. 38]. Hence, this equivalence class is determined by the graph of vertices, or equivalently, by the dual graph of facets. Both are uniquely described by their adjacency matrices. With increasing amount of atoms n, the number of isomers ISO(C n ) grows fast as O(n 9 ), cf. [35]. For instance, there exist a unique C 20 -isomer (we write: |C 20 | = 1) which is dodecahedron, a Platonic solid, but C 60 has already 1812 combinatorial isomers, see [7]. Among all isomers of C 60 , only one (the so-called Buckminster fullerene, an Archimedean truncated icosahedron) has all pentagonal facets being not adjacent. Such isomers are called IPR-fullerenes (from Isolated Pentagon Rule), cf. Figure 1. As illustrated in Figure  2, the number ISO-IPR(C n ) of IPR-isomers also grows asymptotically as O(n 9 ) with increasing number of atoms n [8,32].
The huge variety of possible fullerene isomers with large n (even in the IPR-class) makes the problem of finding molecules with remarkable chemical and physical attributes (including thermodynamic and kinetic stability, permeability, electric conductivity, light diffraction, etc.) extremely difficult. Hence, a need for fast and computationally cheap classification methods of isomers arises. In the literature, there already exist a number of functionals (called chemical descriptors) allowing to find a certain order of isomers, cf. [4,15,20]. These descriptors are of topological, geometric or physical nature. For the practical separation of fullerenes, their potential energetic level is of primary significance. The paper [34] computes the relative energies of all 1812 isomers of C 60 using the density functional theory (DFT) [31] and testing 26 chemical descriptors for their correlation with the energetic ordering of these isomers. The authors identify 7 rules among 26 which they call good stability criteria. These criteria are defined as those able to identify correctly the first two energetically most stable and the three energetically least stable isomers in the correct energetic order such that the correlation coefficient between the ordering according to these rules and the energetic order is at least 0.6. However, the DFT calculations are based on the approximative numerical solution of electronic Schrödinger (linear partial differential) equations requiring from half an hour up to one day of calculation time per isomer (n = 60) on a usual personal computer. Moreover, the convergence of the DFT numerical method is not guaranteed.
Following the famous question by Mark Kac (1966) Can one hear the shape of a drum? [22] we try to "hear" the shape of a fullerene from the spectrum of its adjacency matrix provided that C n -isomers can be mapped bijectively onto their spectra. We also give an example of two different C 44 -isomers with the same spectrum of adjacency matrix of hexagonal facets. We propose a new method of clustering and of classification of C n -isomers for any n ≥ 24, n = 44, based on combinatorial and graph theoretic structure of dual graphs T 6 n of their hexagonal facets. For n = 60 we show in this paper that it yields a good stability criterion with correlation coefficient of 0.9. The spectral analysis of graphs based on adjacency matrices of their vertices has a long standing tradition [3,10,11]. However, we show that for the complete classification and ordering of fullerenes it is sufficient to use • the dual graph T 6 n of hexagons since the positions of 12 pentagons can be reconstructed out of cycles larger than triangles and degrees of vertices in T 6 n (cf. [5]).
• Newton polynomials of the spectrum of the adjacency matrix A 6 n of hexagons up to a certain degree k * since it is well known that they are numerically more stable than the eigenvalues themselves. These polynomials can be computed directly as a trace of A 6 n k , k = 2q, q = 1, . . . , k * /2, whereas the matrix multiplication is computationally less demanding than finding all eigenvalues of a matrix. Hereby, we use a graph theoretical interpretation of Newton polynomials of adjacency matrices of graphs in terms of their cycle numbers.
Our method is computationally very fast requiring O(n 3 log n) operations with a total of 1.15s CPU time on a Intel Core i5-8300H (2.3 GHz) (n = 60). The high correlation of the obtained ordering with the DFT energetic order allows to figure out few energetically stable isomers at a low computational cost. For these isomer candidates, the detailed DFT analysis can be further performed. In order to construct the facet adjacency matrices, a certain enumeration algorithm of all facets is required. Since the spectra of these matrices are invariant with respect to the enumeration of facets, the choice of this algorithm does not matter from the mathematical point of view. For all fullerenes with 24 ≤ n < 380, we use the spiral rule first introduced in [26] (where it is called an orange peel scheme) to enumerate all pentagons and hexagons and generate their adjacency matrices A 5 n and A 6 n , respectively. The first fullerene not obeying the spiral rule is a C 380 -isomer, and the second counterexample is one of over 90 billion C 384 -isomers. All other isomers of C n with n ≤ 450 stick to this rule, cf. [25]. For fullerenes without a facet spiral, a generalized spiral [37] can be used for the one-to-one facet enumeration.
Our spectral approach is illustrated on all isomers of C 60 which are the best studied fullerenes, especially the above mentioned famous Buckminster (soccer ball-like molecule). Such molecular structures are all allotropic forms of carbon [27].

Spectral analysis of C n
For a feasible n a fullerene isomer P ∈ C n is a simple, compact and convex polytope in R 3 with all m := n/2 + 2 facets being either one of 12 pentagons or one of n/2 − 10 hexagons: . . , m. P can be mapped on a two-dimensional graph in a way that edge crossing is avoided and vertex connectivity information is retained. First, one has to choose a facet and rotate P so that this facet is located parallel to the (x, y)-plane at some distance below a fixed projection point q. Next, one draws a line starting in q to each vertex of the polyhedron and extends this line until it crosses the (x, y)-plane. The intersections are the vertices of the new two-dimensional graph, also called Schlegel diagram. Although such a projection is not bijective, it yields a full combinatorial invariant of P . Each of f i can be chosen to be initially parallel to the x − y-plane. Depending on this choice, the resulting graphs can be very different, see [17]. In Figure 3(a) and 3(b), one can see two possible Schlegel diagrams for the Buckminster fullerene (n = 60). The graphs on these Schlegel diagrams are equivalent in the sense that they have the same vertex connectivity. From the definition of a fullerene it immediately follows that the corresponding planar graph is 3-regular. We denote a planar graph of a fullerene by F .  Assume that C n has a facet spiral which is defined as an order of facets such that each facet shares an edge with the previous and next one. This spiral can be presented by a facet spiral sequence. It is a sequence of twelve integers, which determines the position of the twelve pentagons in the facet spiral. Another representation is a sequence of fives and sixes such that the k-th number in the sequence indicates whether the k-th facet in the spiral is a pentagon or a hexagon [1]. We use the first approach in our software [6]. For instance, the facet spiral sequence for Buckminster fullerene represented by its Schlegel diagram in Figure 3(a) is (1,7,9,11,13,15,18,20,22,24,26,32). By C n,i we mean the i th C n -isomer according to the lexicographical order of facet spiral sequences, see also [17, Chapter 2].

Dual facet graphs, adjacency matrices and their spectra
Let G = (V (G), E(G)) = (V, E) be a finite undirected graph with vertex set V and edge set E. Let |V | = m be the number of vertices in V . The adjacency matrix A G = A = (a i,j ) of G is given by In matrix form, A is a symmetric m × m-matrix with zeros on the diagonal and m j=1 a i,j being equal to the valency of the node i: where the inner sum runs over all subgraphs H of G with j nodes and connected components being either edges or cycles, e(H) being the numbers of edges among these components and c(H) being the number of cycles.
b) the Newton polynomial of degree k can be interpreted as the number of all cycles of length k in G.
Here, we call a cycle of length k any closed path (possibly with self-intersections) with k edges from a vertex to itself.
where we put S j = 0, j > m. By [3, Theorem 3.10], it holds b) This interpretation follows immediately from [10, Proposition 1.3.1], since the i th diagonal entry of the k th power of A is the number of walks of length k from vertex i to itself.
For fullerenes, vertex adjacency matrices and their spectra are well-studied, see [1, Section 4.5]. As mentioned above, we generate dual facet graphs T n out of the Schlegel diagrams of C n and consider the spectra of their adjacency matrices. In T n the original facets become vertices, and the original vertices become facets. The edges of the dual graph show adjacency relations between original facets: two nodes of the dual graph are connected by an edge if the corresponding facets of the fullerene are adjacent, i.e. share an edge. In Figure 4(a), one can see the dual graph T 60 of all facets of Buckminster fullerene with the Schlegel diagram in Figure 3(a). Notice that (for the sake of legibility) the facet f 1 is displayed five times in Figure 4(a), whereas it should occur just once. In this paper, we use red, green and white nodes in the images of the dual graphs of C n for pentagons, hexagons and unspecified facets, respectively.
Consider two important induced subgraphs T 5 n and T 6 n of T n . The graph T 5 n illustrates the connectivity between pentagonal facets, i.e. it always contains 12 vertices. For instance, the graph T 5 n of every IPR-isomer consists of 12 disconnected vertices. This being said, it is evident that the number of isomers of C n with the very same graph T 5 n increases rapidly with increasing n, cf. Figure 2 for the IPR-case. Hence, considering the graph T 5 n does not yield an invariant for all C n -isomers. In order to characterize all isomers we need the graph T 6 n showing the connectivity of all m 6 hexagonal facets of a C n -isomer. As an example the graph T 6 60 of the Buckminster fullerene is shown in Figure 4(b). It turns out that T 6 n completely characterizes the graph T n . We denote by A n , A 5 n , A 6 n the adjacency matrix of T n , T 5 n and T 6 n , respectively. Remark 1. a) For k > m a formula similar to (1) can be derived, by which it follows that traces of the k th power of A with k > m are linear combinations of traces of smaller powers. This can be explained by the fact that subgraphs have at most as much vertices as the whole graph.
b) An alternative approach is to insert formula (2) into itself, which yields a representation of Newton polynomials as a polynom with several unknowns being Newton polynomials of lower degree.
c) The condition that every vertex in a fullerene has valency three corresponds to the fact that the dual graph T n consists of triangles only. However, the subgraphs T 5 n and T 6 n may also have larger cycles. d) No 4-cycles exist neither in T n nor in T 5 n and T 6 n , see [13,Theorem 4.15 (1)]. e) The problem of description of all simple cycles (i.e., closed loops without self-intersections) of pentagonal or hexagonal facets of length k is crucial to combinatorial classification of fullerenes. Denote by N 0 the set of natural numbers and zero. Obviously, it holds N (A 6 n , k) ∈ N 0 for all k ∈ N. For the dual graphs T 5 n and T 6 n we construct their adjacency matrices A 5 n and A 6 n . Since the number of unit entries in each line does not exceed 5 for A 5 n or 6 for A 6 n it holds by Gershgorin's theorem that In the sequel, let us concentrate on the properties of σ(A 6 n ). The classical Frobenius-Perron theory applied to graph spectra [10, Proposition 3.1.1] yields a more accurate estimate for the largest eigenvalue of A 6 n which is positive and of multiplicity one if T 6 n is connected.
Lemma 2. Let G be a graph with m vertices and valencies κ 1 , . . . , κ m , and H be an induced subgraph of G.
Let λ max (G) and λ max (H) be the largest eigenvalues of the adjacency matrix A G and A H . a) If G is connected, then κ min ≤κ ≤ λ max (G) ≤ κ max , withκ being its mean valency and κ max the maximum valency of the vertices in G. In particular, Notice that for all connected irregular graphs T n , T 5 n and T 6 n we haveκ > 1 and κ max ≤ 6. For nonconnected or irregular graphs we get only 0 ≤κ ≤ κ max . For instance, the dual graph T 6 60 of the Buckminster fullerene (see Fig. 4(b)) is connected and regular with κ = 3, thus the largest eigenvalue of A 60,6 for it is equal to 3.
In general, the value θ := κ max −κ ≥ 0 is a measure of the asymmetry of the dual facet graph which we call the asymmetry coefficient.
is the indicator function of B and a max is the multiplicity of the eigenvalue λ max . Indeed, by [3, Theorem 6.3], it holds −λ max ∈ σ(A 6 n ) iff our dual graph T 6 n is bipartite, i.e., it has no cycles of all odd lengths, cf. [3,Corollary 3.12]. In this case, −λ max has necessarily multiplicity one, see [3,Lemma 3.13]. Since for large n ≥ 40 the dual graph T 6 n of hexagons of any isomer of C n contains either a 3-or a 5-cycle, it holds −λ max ∈ σ(A 6 n ), and the behavior of the the whole Newton polynomial tr A 6 n k for large powers k is dominated by λ k max .
Let |A| be the cardinality of a finite set A.
Lemma 3. For all feasible n, consider a C n -isomer with graphs T n , T 5 n and T 6 n . Then it holds Proof. Since the graph T n of a C n -isomer is 3-regular, it has 3n 2 edges. An edge (v, w) ∈ E(T n ) is either an edge between two pentagons, (v, w) ∈ E(T 5 n ), or between two hexagons, (v, w) ∈ E(T 6 n ), or between a pentagon and a hexagon, (v, w) ∈ E(T n \ T 5 n ∪ T 6 n ). Since the number of pentagons is always 12, there are 60 − 2|E(T 5 n )| edges between pentagons and thus which finishes the proof.

Cospectral Isomers
It is easy to prove that isomorphic graphs are cospectral. The inverse statement is in general not true (cf. e.g. a counterexample in [36]). The natural question arises: Which graphs are determined by their spectrum? [36]. Some specific graphs like paths, complete graphs, regular complete bipartite graphs, cycles and their complements yield a positive answer to this question. In what follows, we provide new examples of non-isomorphic cospectral graphs coming from the world of fullerenes.
We derive some theoretical results for sets of non-cospectral graphs and apply them to fullerene isomers. To begin with, recall the following Since For all even 24 ≤ n ≤ 150, we checked the existence of cospectral pairs of C n -isomers within T n , T 5 n and T 6 n graphs. Using Lemma 4, we applied MATLAB functions eig(), charpoly() and trace() with double precision, cf. [28]. In cases where two isomers seem to be cospectral with respect to T n or T 6 n , we increase the precision by using vpi format and MATLAB symbolic toolbox [28].
For n < 32 no pair of cospectral isomers with respect to T n , T 5 n and T 6 n can be found. Our results for 32 ≤ n ≤ 60 are listed in Table 1. This table can be extended to n > 60 with all zeroes in its first and third rows, and positive integers in its second row.
Considering the whole dual graph T n one finds only one pair of cospectral isomers with n = 44, compare Figures 5(a) and 5(b). To explain the difference within this pair, we need the following Definition 2. Let G and H be two graphs. a) A fragment F of G is a connected induced subgraph of G. In particular, for fullerenes we call fragments of T 5 n and T 6 n pentagon-fragments and hexagon-fragments of T n , respectively. b) Two non-isomorphic graphs G and H are called fragment flipped if two isomorphic fragments F G in G and F H in H exist such that the remaining graphs G \ F and H \ F are isomorphic.
As one can see in Figure 5(a) and 5(b), the two cospectral C 44 -isomers are fragment flipped. However, in general such a flip does not preserve the spectrum of a graph. To illustrate this, Figure 5(c) and 5(d) contains a non-cospectral fragment flipped pair of C 60 -isomers.
For n ∈ {32, 36, 40, 52} a pair of distinct isomers exists with the same spectrum σ T 6 n , since their graphs T 6 n are isomorphic.  Table 1: Number of non-unique spectra of graphs T n , T 5 n and T 6 n of C n -isomers.
For a fixed 40 ≤ n ≤ 150 at least three and at most 298 non-unique spectra σ(T 5 n ) exist. Since the amount of different arrangements of 12 pentagons is limited and the number of isomers grows rapidly with increasing n, there must exist isomers with same subgraphs T 5 n . We discuss the number and different configurations of pentagon-fragments in [5] in more detail.
The above empirical findings lead to the following Conjecture 1. a) For all feasible n = 44 two C n -isomers are isomorphic iff they are cospectral with respect to T n .
b) For all feasible n ≥ 54 two C n -isomers are isomorphic iff they are cospectral with respect to T 6 n .
c) For any feasible n at least one of the spectra σ(T n ), σ(T 5 n ), σ(T 6 n ) is unique for all C n -isomers.
Proof. Ordering all absolute spectra |σ(G)| lexicographically yields a unique order on the set Γ. Then all graphs G from Γ can be distinguished either by |σ [1,∞) (G)| or by |σ [0,1) (G)|. In the first case, let us consider all G ∈ Γ with distinct |σ [1,∞) (G)|. The sum λ∈|σ [0,1) (G)| |λ| k converges to 0 for k → ∞, so its influence on the Newton polynomials N (A G , k) for k large enough can be neglected. Since |σ [1,∞) (G)| is unique for all considered graphs G there must exist a degree k * 1 such that the values of N (A G , k) for all even k ≥ k * 1 are distinct. Now if some graphs G ∈ Γ have the identical part of absolute spectrum |σ [1,∞) (G)|, they must differ within the part |σ [0,1) (G)|. Hence, there must always exist a degree k * 2 ≤ k * 1 such that the sum λ∈|σ [0,1) (G)| |λ| k * 2 distinguishes these graphs. Lemma 4 c) and Remark 1 a) yield the number of vertices in the graph G as an upper bound for k * 1 and k * 2 .
Corollary 1. Assuming Conjecture 1 to be true, all C n -isomers can be uniquely characterized by at most two Newton polynomials N (A G , k * 1 ) and N (A G , k * 2 ), k * 1 ≤ k * 2 even, with respect to at least one of the graphs G ∈ {T n , T 5 n , T 6 n }.
Proof. This follows from Theorem 1 and Lemma 2, since λ max (G) > 1 for at least one of the graphs G ∈ {T n T 5 n , T 6 n }.
Remark 3. (i) Note that the set of C n -isomers with distinct largest eigenvalues λ max (G) for some G ∈ {T n T 5 n , T 6 n } has distinct sets |σ [1,6] (G)|.
(ii) It is important to stress that we consider even degrees k ≤ m only. In general, values of Newton polynomials N (A G , k) with odd degrees k do not distinguish between graphs G. To illustrate this point, consider the hexagon graphs of two specific isomers of C 28 . The graph T 6 28,1 of the first isomer consists of four isolated hexagons and the other graph T 6 28,2 of two pairs of two adjacent hexagons. One gets the following spectra: σ(T 6 28,1 ) = {0, 0, 0, 0}, σ(T 6 28,2 ) = {−1, −1, 1, 1}.
It follows directly that Newton polynomials of any odd degree do not distinguish between C 28,1 and C 28,2 , but N (A 6 28,1 , k) = N (A 6 28,2 , k) for every even k ≥ 2.

Spectral classification and stability prediction
In this section, we decompose the family of all C n -isomers into subsets, which we call clusters, using the Newton polynomials N (A G , k), k = 2, 4, 6, . . . , k * with adjacency matrix A G . This approach can be applied to any graph G ∈ {T n , T 5 n , T 6 n } with no cospectral isomers. Definition 3. a) For a feasible n and a given even integer k, we call a set of isomers with same value N (A G , k) of Newton polynomial of degree k a cluster of C n . For a fixed n, the family of all these clusters is called a clusterization of C n . A clusterization such that its every cluster has exactly one element is called a complete clusterization. b) We define k * single as the minimal degree k of Newton polynomials which is needed for a complete clusterization.
As we have seen in Lemma 1 b), the Newton polynomial of degree two and three is equal to twice the number of edges and six times the number of triangles in the considered graph. The interpretation of N (A n , 4) is a bit more complex. We have to count all possible cycles of length four. In Figure 6, the idea of calculation of N (A 60 , 4) is illustrated on Buckminster fullerene C 60,1812 . There, all five possible cycles of length four together with their frequencies are listed for one pentagon in C 60,1812 .  For the graph T n , one can show that the sum of these frequencies over all vertices in T n only depends on n, cf. [5]. So, N (A n , k) for k ≤ 4 can be neglected for the clusterization of C n on the basis of T n , since these Newton polynomials have the same value for all C n -isomers.
Nevertheless, for T 6 n we have to consider every Newton polynomial of even degree, since in T 6 n the number of vertices is fixed, but neither the number of edges between them nor the number of triangles in T 6 n is determined by n.
In the following section we use T 6 60 in order to get a complete clusterization of C 60 as an example.

Clusterization of C 60 using Newton polynomials
Recall that by Euler's formula each C 60 -isomer has 90 edges and 32 facets with 12 pentagons and 20 hexagons among them. The Newton polynomial N (A 6 60 , 2) can take on 18 distinct values 60 = t 1 < . . . < t 18 = 100. For values t 1 = 60, t 2 = 64, t 16 = 92, t 17 = 96 and t 18 = 100, there exists exactly one C 60 -isomer with N (A 6 60 , 2) = t i , i ∈ {1, 2, 16, 17, 18}. These isomers are C 60,1812 , C 60,1809 , C 60,2 , C 60,3 , C 60,1 , respectively. Their Schlegel diagrams and dual hexagonal graphs T 6 60 are shown in Figure 7. Moreover, our numerical results show that for any degree k ≥ 2 with respect to T 6 60 these isomers form a cluster with one single element. Ordering these five isomers according to N (A 6 60 , k) does not change with increasing degree k ≥ 2. In addition, Newton polynomials of all other C 60 -isomers are bounded by N (A 6 60,1809 , k) and N (A 6 60,2 , k), i.e. N (A 6 60,i , k) ∈ N (A 6 60,1809 , k), N (A 6 60,2 , k) for all i / ∈ {1, 2, 3, 1809, 1812} and all k ≥ 2. We checked that for any pair of two C 60 -isomers with distinct Newton polynomials of degree k the Newton polynomials N (A 6 60 ,k) with k ≤k ≤ 100 are distinct as well. So, it holds k * 1 = k * 2 = k * single = 12, where k * 1 and k * 2 are from Theorem 1. One can observe that the number of distinct Newton polynomials and so the number of clusters with one element is monotone growing with k. Numbers of clusters and clusters with one element for all even 2 ≤ k ≤ 100 are listed in Table 2.  Table 2: Number of clusters with respect to Newton polynomial N (A 6 60 , k) for even 2 ≤ k ≤ 100 We applied the above clusterization scheme to C n , 28 ≤ n ≤ 150 and plotted n against k * single in Figure 8(a). Recall that a (pessimistic) upper bound for k * single is the number of vertices in T 6 n , i.e. k * single ≤ m 6 = n 2 − 10 due to Lemma 4. However, the good news is that the actual growth rate of k * single is logarithmic with n. Using MATLAB curve fitting toolbox [28] we get Figure 8: Minimal degree k * single (a) and k * pair (b) needed for the complete clusterization of C n -isomers with 28 ≤ n ≤ 150 Next, we use pairs of Newton polynomials N (A 6 n , k 1 ), N (A 6 n , k 2 ) with k 1 < k 2 ≤ k * single in order to cluster all C 60 -isomers. By considering additionally a second Newton polynomial of lower degree, we hope to decrease the needed degree to get a complete clusterization. We define k * pair as the minimal k 2 such that a complete clusterization is given. Indeed, this approach decreases the needed degree significantly (compare both plots in Figure 8), and therefore, reduces computational costs. For C 60 and T 6 60 the following four tuples with k 2 < k * single of degrees of Newton polynomials lead to a full classification: (4,10), (6,10), (8,10) .
We get k * pair = 8. Next we plotted all values for k * pair against n and assumed a logarithmic function as for k * single . Using MATLAB curve fitting toolbox we get the following approximation k * pair (n) ≈ −24.83 + 10.29 log (0.334n + 9.989) with a coefficient of determination of R 2 = 0.8277. A third hierarchical approach in order to decrease the needed degree uses a vector with all Newton polynomials up to degree k ≤ k * single . Analogously to the first two approaches, we define k * hierarchical as the minimal k which yields a complete clusterization. This approach decreases e.g. the degree for n = 78 from k * pair = 12 to k * hierarchical = 10. So, this approach does not change the needed degree significantly. Nevertheless, we performed the same interpolation using MATLAB curve fitting Toolbox and got k * hierarchical (n) ≈ −97.05 + 19.83 log (1.466n + 125.5) .
with R 2 = 0.8818. In [5] we discuss another upper bound for k * single (sharper than m 6 ) using the generalized Stone-Wales operation introduced in [2]. In the paper, we give a combinatorial interpretation of k * and derive an equation system with Newton polynomials which determines whether a fullerene with given n can be constructed.

Relative energy
A fullerene isomer, which can be chemically separated with a significant mass quantity and uniquely characterized, is called stable. In order to decide which C 60 -isomer can be stable the relative energy of all of them was calculated with high-accuracy quantum chemistry methods and discussed in [34]. Here, relative means compared with the Buckminster fullerene C 60,1812 which has the lowest DFT-energy at the PW6B95-D3 ATM /def2-QZVP level (cf. [34]), i.e. Buckminster fullerene has a relative energy of 0. In the sequel, we say that an isomer C n,i is energetically more stable than C n,j , i = j, if C n,i has a smaller energy than C n,j .
According to [34] the most stable isomer is C 60,1812 and the second stable one is C 60,1809 . At the other end of the ranking the three least stable ones are C 60,2 , C 60,3 and C 60,1 . It has been assumed for a long time that an isomer is the more stable the less adjacent pentagon it has. Indeed, calculations of [34] allow the conclusion that each pair of two adjacent pentagons leads to a increase in the relative energy of an isomer of about 20 to 25 kcal mol −1 . The amount of such pentagon pairs can be described with Fowler-Manolopoulos pentagon indices p i := #{pentagons which are adjacent to i other pentagons}, such that the sum of p 1 up to p 5 is equal to 12 for every fullerene, cf. [17]. Based on these values, the pentagon signature P 1 = 1/2 5 i=1 ip i can be calculated, which quantifies the amount of connected pentagons. Clustering all C 60 -isomers according to the pentagon signature, five isomers stand out, namely C 60,1812 (P 1 = 0), C 60,1809 (P 1 = 2), C 60,2 (P 1 = 16), C 60,3 (P 1 = 18) and C 60,1 (P 1 = 20). The signature can be easily read from Figures 7(a)-7(e). Pentagons signatures of the remaining isomers lie between 2 and 16. For each of these values, at least two isomers exist with the same pentagon signature. By Lemma 1 b) and 3 we get the following Proposition 1. For any C n -isomer it holds In Table 3, C 60 -isomers are listed in the same order given by their relative energy, by their pentagon signature and their Newton polynomial of degree 2.
Indeed, one gets more information about a fullerene structure looking on hexagons than on pentagons. This becomes clear looking at fullerenes with large n. For example, C 80 has 7 IPR-isomers, so their pentagon structure and pentagon signature are the same. Nevertheless, only two of them have been produced in pure form, although DFT calculations have been done for all of them, see [23]. As a result of [23], only two IPR-isomers can be claimed stable. Hence, all descriptors based on the pentagon structure do not properly predict stability.
In [34] a good stability criterion is defined as the one which can identify C 60,1812 and C 60,1809 as the most stable and C 60,2 , C 60,3 and C 60,1 as the least stable isomers in the correct energetic order. Additionally, the Pearson coefficient ρ of linear correlation between the relative energies of all C 60 -isomers and their criterion values should be larger than 0.6. Finally, the slope and the Pearson correlation coefficient in the linear regression of relative energy vs. the criterion for C 60 -isomers with P 1 ∈ {4, . . . , 14} should have the same sign.
As we have seen in Table 3, Newton polynomials yield the correct order of the most and least stable C 60 -isomers. Next, we perform a linear regression (using MATLAB curve fitting Toolbox) of N (A 6 60 , k) vs. relative energies of all C 60 isomers for all even 4 ≤ k ≤ 12. For the case k = 2 Newton polynomial is equivalent to the 1st moment hexagon Signature H 1 , which is listed in [34, Table 3] as a good stability criterion. Hence, N (A 6 60 , 2) is a good stability criterion as well and can be neglected in further considerations. Table 4 shows that Pearson correlation coefficient ρ is much higher than 0.6 for all considered k. For degrees k = 8, 10, 12, Newton polynomials N (A 6 60 , k) get very large, and therefore we took a logarithmic scale. But even with linear scale, one gets correlation coefficients larger than 0.6 in these three cases, compare Table 7.
Next we divided all isomers of C 60 into 18 subsets G i = {P ∈ C 60 | P 1 (P ) = i} according to their pentagon signature i ∈ {0, 2, 3, . . . , 16, 18, 20}. Then we performed a linear regression of Newton polynomials of different degrees vs. relative energies of isomers in G i for every i / ∈ {0, 2, 3, 15, 16, 18, 20} as it is required in [34]. Our results are listed in Appendix, Table 8. For k ∈ {4, 6} and i ∈ {4, 5} we get Pearson correlation coefficients and slopes with a negative sign, unlike for all other combinations of k and i. This can be explained by the fact that G 4 and G 5 do not contain many isomers. More precisely, |G 4 | = 17, and |G 5 | = 86 holds. So, neglecting these two cases would yield that N (A 6 60 , k) with k = 4, 6 is a good stability criterion. For k = 10, 12 one gets positive slopes and Pearson correlation coefficients in all cases, and therefore Newton polynomials of degree 10 and 12, in particular of degree k * single , entirely fulfil all conditions of a good stability criterion.  Table 4: Pearson correlation coefficient ρ and the slope of linear regression between Newton polynomials and the relative energy of all C 60 -isomers given in [34].
To check whether Newton polynomials can distinguish between IPR-isomers, i.e. yield their energetically correct order, we computed Newton polynomials of all 31924 C 80 -isomers. Within the whole set of C 80 the seven IPR-isomers have the smallest Newton polynomials. But ordering the set of IPR-isomers according to Newton polynomials leads to the observation that the most stable IPR-isomers have the greatest Newton polynomials. These seven isomers are listed in Table 5. So, it seems that with increasing n Newton polynomials N (A 6 n , k) for even k ≥ 2 remain a good stability criterion.

Asymmetry coefficients of isomers of C 60
The Fowler asymmetry parameter is claimed to be a good stability criterion [34]. Check whether the asymmetry coefficient θ defined in Section 2 is a good stability criterion as well. The asymmetry coefficients of isomers shown in Figure 7 are θ 1 = 1, θ 2 = 1.4, θ 3 = 1.2, θ 1809 = 0.8 and θ 1812 = 0. The histogram of θ of all C 60 -isomers is shown in Figure 9. It turns out that the asymmetry coefficient θ is not a good stability criterion since it does not preserve the energetic order required in [34]. Thus, the asymmetry coefficient of eight isomers is equal 2.4, which is the largest value. These isomers are C 60,1334 , C 60,1554 , C 60,1676 , C 60,1740 , C 60,1741 , C 60,1742 , C 60,1761 and C 60,1784 . Figure 10 shows C 60,1784 , which has six hexagons with valency two, twelve hexagons with valency four and two hexagons with valency six. So, the mean valency is 3.6, the maximal valency is 6 and the resulting asymmetry coefficient equals 2.4. In Table 6, asymmetry coefficients and energetic order numbers (1 = most stable, 1812 = least stable) [34] are listed for some outstanding isomers of C 60 , i.e. the ones with smallest and largest θ as well as those shown in Figure 7. One can see that the three least stable as well as the two most stable isomers have the asymmetry coefficient less than 2.4. In addition, the relative stability of the eight isomers with θ = 2.4 varies from 313 to 1122. This leads to the conclusion that θ is not a good stability predictor.  Table 6: Comparison of the asymmetry coefficient θ and energetic stability of some C 60 -isomers

Conclusion
We present an easy to compute functional of spectra of the graphs T n and T 6 n which classifies all C n -isomers. Thereby we focus on the structure of the dual graph of hexagonal facets of C n and its adjacency matrix A 6 n . The spectra of the adjacency matrices are characteristic to combinatorial isomers described above. It becomes apparent that the Newton polynomial of degree 2, 10, 12(= k * single ) of T 6 n appears to be a good stability criterion. So, Newton polynomials of T 6 60 can be added to the list presented in [34, Table 3] as indices, which, depending on the degree, fulfil the criteria partly or entirely. We show that Newton polynomials generalize the Pentagon signature and better describe the fullerene structure. The interpretation of these Newton polynomials is very easy for k ≤ 3, but gets demanding with increasing k.    Table 8: Linear regression of an independent variable N (A 6 60 , k), 4 ≤ k ≤ 12 even, vs. relative energy over the subsets G i of C 60 -isomers