The Interplay between Diradical Character and Stability in Organic Molecules

: The number of scientiﬁc papers on the unique properties and the potential for various applications of compounds with a diradical character is growing constantly. The diradical character enhances and even engenders certain desired optical properties and its modulation is a modern molecular design strategy. Nowadays, molecules with a non-zero diradical character are regarded as promising materials for new-generation and highly efﬁcient solar cells and photonics devices. What is the price, however, of the unique properties of open-shell compounds? Alongside all the beneﬁts, the diradical character is usually associated with low stability and high reactivity—unwanted molecular qualities for practical purposes. Thus, from a fundamental and applied point of view, it is important to investigate the correlation between the diradical character and laboratory stability, which is the goal of the present paper. Here, we report a combined quantum–chemical study (conceptual DFT and spin-projected HF theory) and multivariate analysis of the diradical character of a series of o and p -quinomethides, for the stability of which experimental data are available. Our results reveal that a compromise between the diradical character and laboratory stability of a molecule is feasible and that the relationship between these two quantities can be understood in the framework of Clar’s sextet theory.


Introduction
The diradical character (DRC) determines the chemical reactivity and stability of molecules and produces their unique optical, magnetic, and electronic properties [1]. It varies between zero (for closed-shell systems) and one (for diradicals) and is related to the fundamental understanding of the chemical bond. The diradical character was originally defined as twice the weight of the doubly excited configuration in the singlet ground state wavefunction within the multiconfigurational self-consistent field theory and is not a directly observable feature experimentally [2]. In 2010, the connection between the DRC and spectroscopically measured quantities was demonstrated, which allowed the experimental estimation of the DRC [3]. To the best of our knowledge, however, in the literature, there are few examples of DRC evaluated by spectroscopic measurements.
Recently, it was demonstrated that the tunability of the DRC represents a fruitful computational approach for the discovery of new organic materials [4][5][6][7][8]. For example, molecules with a moderate (0.5) diradical character, i.e., diradicaloids, possess unusually high values of the second hyperpolarizability and are attractive for application in organic photonics [9][10][11][12]. It was also shown that compounds with a low-to-intermediate (0.1-0.5) diradical character are potential singlet fission materials [13][14][15][16]. Singlet fission is a phenomenon in which one absorbed photon is converted into two excitons. When implemented in a solar cell, it can double the number of charge carriers [17]. Therefore, Symmetry 2021, 13, 1448 2 of 11 molecules able to transform one photon into two excitons are regarded as next-generation, highly efficient organic photovoltaic materials [18][19][20][21]. From a thermodynamics point of view, in singlet fission chromophores, twice the excitation energy to the first excited triplet state, E(T 1 ), should be less or equal to the excitation energy to the first excited singlet state, E(S 1 ): 2E(T 1 ) − E(S 1 ) ≤ 0. It was shown that the thermal losses in singlet fission are proportional to l2E(T 1 ) − E(S 1 )l and, hence, zero or slightly smaller differences are required from a practical point of view [4]. Since the first triplet excitation energy is also inversely proportional to the DRC, it is not surprising that molecules with a DRC between 0.1 and 0.5 provide optimal 2E(T 1 ) − E(S 1 ) values for real photovoltaics application.
Alongside attractive properties, the diradical character is usually an implication for high molecular reactivity and instability [22]. Here we arrive at an intrinsic paradox: the presence of the diradical character, which makes a substance prospective for photovoltaics and non-linear optics, may also disqualify the material as inappropriate for application. A deeper analysis of this contradiction reveals that it expresses the balance between the properties of the ground and the excited states.
Molecular stability is usually understood in two aspects: thermodynamic and kinetic. The thermodynamic stability is determined solely by the molecular structure and is related to the Gibbs energy, which is governed by the total electron energy of the molecule [23]. Thus, the diradical character depends on the HOMO-LUMO gap (where HOMO is the highest occupied molecular orbital and LUMO is the lowest unoccupied molecular orbital) and, as an electronic characteristic, is associated with the thermodynamic stability of the compound. The kinetic stability is also structure-dependent and is defined by the chemical reactivity of the molecules and the rates of the reactions in which they participate [23]. The chemical reactivity depends on the ionization potential, the electron affinity, the charge and spin distribution (electron density) in a molecule and, hence, is correlated with its diradical character as well. A kinetically unstable (highly reactive) molecule may have a fairly long lifetime if kept away from agents for which it has affinity (e.g., in a vacuum). A thermodynamically unstable molecule, however, even if completely isolated, has a short lifetime because, to stabilize, it will undergo structural changes (e.g., Jahn-Teller effect) or it may never become synthetically feasible.
To date, the relationship between reactivity, stability, and the diradical character remains an intriguing topic from a fundamental perspective [22]. The present theoretical study aims to contribute to the deeper understanding of the problem and to explore whether there are stability limitations to the application of molecules with a non-zero DRC. For this purpose, we combine quantum-chemical and chemometrics approaches to treat a family of quinone methides, which differ in laboratory stability and reactivity. To the best of our knowledge, this is the first theoretical study to focus on the intimate relationship between the diradical character and laboratory stability of molecules.

Molecular Set
The study is based on an old, but fruitful for the purpose, paper of Musil et al. who tried to explain the laboratory stability of quinone methides applying Hückel theory, i.e., with respect to aromaticity [24]. A valuable contribution of the authors is the summary of the literature data on the laboratory stability of a large number of structurally similar but markedly differing in reactivity conjugated compounds ( Figure 1). Presented are also the mechanisms of the reactions in which the quinone methides participate, causing their laboratory instability. The occurrence of an exocyclic methide group most often initiates reactions with nucleophiles (water) or with a radical mechanism (dimerization). The authors define as stable those quinone methides that can be isolated and have infinite lifetime, allowing detailed examination of their properties. Moderate stability applies for quinone methides that can exist only in solution or have a limited lifetime (minutes or hours) but long enough to allow spectroscopic validation and investigation. As unstable are classified the highly reactive quinone methides with too short lifetimes, which makes their experimental characterization impossible. Some analogous compounds differ only in the type of substituents and therefore have two versions labelled as 'a' for the hydrogen substituent and 'b' for the phenyl-group substituent. Those that have no analogues bear no label and have only hydrogen substituents in the exocyclic methide group.
Symmetry 2021, 13, x FOR PEER REVIEW 3 of 12 time, allowing detailed examination of their properties. Moderate stability applies for quinone methides that can exist only in solution or have a limited lifetime (minutes or hours) but long enough to allow spectroscopic validation and investigation. As unstable are classified the highly reactive quinone methides with too short lifetimes, which makes their experimental characterization impossible. Some analogous compounds differ only in the type of substituents and therefore have two versions labelled as 'a' for the hydrogen substituent and 'b' for the phenyl-group substituent. Those that have no analogues bear no label and have only hydrogen substituents in the exocyclic methide group.

Calculation of the Diradical Character
There are two general approaches in the theoretical prediction of DRC: the use of multiconfigurational methods such as RASSCF/CASSCF [25] or the application of the spin-projected unrestricted Hartree-Fock (PUHF) method [2] with natural orbital analysis [26]. The multiconfigurational methods are highly accurate for open-shell systems and they can correctly predict the DRC of a singlet fission chromophore, as well as all properties that determine the application potential of a molecule. However, the multiconfigurational methods are heavily resource-demanding and usually limited to small or medium molecules. The PUHF method allows finding the broken symmetry singlet solution and, in combination with natural orbital analysis, can be used for DRC calculations even for very large molecules at a negligible computational cost. PUHF is an approximated method for the treatment of open-shell systems and, while it accurately predicts the DRC, it cannot estimate a variety of other molecular properties. In the framework of the PUHF theory, the diradical character (y0, i = 0) depends on the occupation numbers of HONO (highest occupied natural orbital) and LUNO (lowest unoccupied natural orbital):

Calculation of the Diradical Character
There are two general approaches in the theoretical prediction of DRC: the use of multiconfigurational methods such as RASSCF/CASSCF [25] or the application of the spin-projected unrestricted Hartree-Fock (PUHF) method [2] with natural orbital analysis [26]. The multiconfigurational methods are highly accurate for open-shell systems and they can correctly predict the DRC of a singlet fission chromophore, as well as all properties that determine the application potential of a molecule. However, the multiconfigurational methods are heavily resource-demanding and usually limited to small or medium molecules. The PUHF method allows finding the broken symmetry singlet solution and, in combination with natural orbital analysis, can be used for DRC calculations even for very large molecules at a negligible computational cost. PUHF is an approximated method for the treatment of open-shell systems and, while it accurately predicts the DRC, it cannot estimate a variety of other molecular properties. In the framework of the PUHF theory, the diradical character (y 0 , i = 0) depends on the occupation numbers of HONO (highest occupied natural orbital) and LUNO (lowest unoccupied natural orbital): where T i is the overlap integral of the respective orbital pair, which is calculated using the occupation numbers (n) of the natural orbitals of the UHF wavefunction:

Quantum-Chemical Protocol
Our computational protocol encompasses geometry optimization of the molecules with B3LYP/6-31G*, followed by frequency analysis [27][28][29][30]. Next, the optimized structures are subject to single point "broken symmetry" calculations with the PUHF/6-31G** method in combination with natural orbital population analysis, from which the diradical character is quantified. All calculations were performed with Gaussian 09 [31]. The chemical reactivity, the global and local reactivity indices are calculated within the conceptual density functional theory (CDFT) [32]. CDFT is a popular and well-proven approach for the extraction of chemically relevant principles from DFT, in particular for interpretation of the chemical reactivity of the molecules. All CDFT data reported in the manuscript were obtained using the program Multiwfn 3.7 [33]. In order to calculate the CDFT quantities for a given molecule, the Multiwfn 3.7 program requires as input three separate wavefunctions with molecular charges of 0, +1, and −1. These wavefunctions were obtained from B3LYP/6-31G* single-point calculations on the optimized geometry with the Gaussian 09 software. The computed values of the local and global reactivity indices and their definitions are collected in Table 1. For more details on the approximations used for the calculation of all CDFT quantities, the readers are referred to the manual of the Multiwfn 3.7 program and all related references therein. Table 1. Global and local reactivity indices computed within the CDFT framework [33]. N is the number of electrons in a molecule.

Global (B3LYP/6-31G*)
Local (B3LYP/6-31G*) with Hirschfeld Charges Vertical electron affinity Largest atomic Fukui index for radical attack 2η Electrophilicity In addition to the computed global and local reactivity indices, we adopted all quantities obtained by Musil et al. with Hückel theory, namely, REPE (resonance energy per π-electron), ∆E i , ∆E ij (reactivity indices), and P (stability index) [24]. The stability index is a combined quantity that we derived to explain the laboratory stability of quinone methides and is calculated as follows: Based on the stability index, the authors [24] have succeeded to predict and synthesize moderately stable (IIb, Vb and VIb) and unstable (VIa) quinone methides.
Complementing Hückel theory, another iconic rule is used for interpretation of the results, namely, Clar's rule, stating that: the Kekulé resonance structure with the largest number of disjoint aromatic π-sextets is the most important for characterization of the properties of polycyclic compounds; the aromatic π-sextets are defined as six π-electrons localized in a single benzene-like ring separated from adjacent rings by formal C-C single bonds [34].

Multivariate Analysis
The combination of multivariate statistical methods with atomistic simulation and modeling has become a useful tool in the field of structure-property classification and prediction models [35][36][37]. The classification model used in this study is based on the grid of property descriptors described in the previous section. The hierarchical cluster analysis (HCA) and a non-hierarchical (K-means) method were used as exploratory data analysis methods to reveal the main separation and discriminating factors for the dataset.
As an unsupervised clustering method, HCA is employed to gain a better insight into data structures; the formed groups (clusters) are based on a similarity measure (i.e., Euclidean distance) without considering the class-membership function [38,39]. The graphical representation of the clustering results can be displayed in a tree-like structure, the dendrogram, which can be used to probe and recognize the similarity patterns and clusters within objects and variables in a feature space.
To conduct unsupervised learning, non-hierarchical clustering is applied. K-means clustering is one of the common unsupervised machine learning algorithms, applied for a plethora of problems. The K-means algorithm classifies n observations into k clusters grouped around centroids by assigning every data point from the input matrix to the nearest cluster, keeping the distance to the centroids as short as possible. The number k denotes the centroids, and the centroid is the real or imaginary location indicating the center of the cluster. Combined sequentially, the two approaches can efficiently use the local descriptor environment generated by quantum chemistry calculation and can produce a satisfying and chemically accurate classification model revealing patterns of similarity between the objects of study. All chemometric calculations were performed with the software package STATISTICA 8.0 [40].

Results and Discussion
The results of our study on quinone methides are presented in Figure 2, while Figure 3 summarizes the results reported by Musil et al. [24]. All numerical data are presented in the Supplementary Materials (SM), Table S1. Figure 2 clearly illustrates that all the quinone methides considered have a non-zero diradical character. DRC varies in the series in the range 0.07-0.51 and, consequently, with respect to the requirements for singlet fission chromophores, these compounds are potential photovoltaics. The results for all stable quinone methides isolated in the crystal phase are coloured in green and they feature DRC values in the interval of 0.07-0.20. For comparison, the respective values for the unstable species vary within the range of 0.11-0.40. Consequently, the laboratory stability of quinone methides cannot be explained solely with the DRC values. Particularly illustrative is compound VIa, featuring one of the highest DRC values in the series, notwithstanding exhibiting moderate laboratory stability. Another interesting observation is that the functionalization of the exocyclic carbon with two phenyl substituents results, in most cases, in augmented DRC (by 0.08 on average) concurrent with enhanced kinetic stability.
Among the global reactivity indices, the global hardness mirrors the DRC values, while the softness changes synchronously with the DRC. Considering the definitions of hardness and softness (Table 1), this is an expected result: molecules with high softness/low hardness have a narrow HOMO-LUMO gap, hence, a pronounced diradical character.
It is important to point out the parallel between the diradical character and the electron affinity. The compounds with strong electron affinity, i.e., with low LUMOs, are those with a large DRC. In the language of chemistry, strong electron affinity signifies an easy reduction. This relationship is not unexpected, since, in some alternative approaches for DRC assessment, the DRC is defined as the occupation number of LUNO, which increases with the decrease in LUNO energy. Correspondingly, the low ionization potential indicates a high HOMO, and this, combined with a low LUMO, results in a narrow HOMO-LUMO gap and a substantial DRC. In most cases, the low ionization potential is related to a large DRC, but the ionization potential alone cannot explain the variations in the DRC. The electronegativity essentially follows the profile of the ionization potential, while the electrophilicity matches softness.    It is important to point out the parallel between the diradical character and the electron affinity. The compounds with strong electron affinity, i.e., with low LUMOs, are those with a large DRC. In the language of chemistry, strong electron affinity signifies an easy reduction. This relationship is not unexpected, since, in some alternative approaches for DRC assessment, the DRC is defined as the occupation number of LUNO, which increases with the decrease in LUNO energy. Correspondingly, the low ionization potential indicates a high HOMO, and this, combined with a low LUMO, results in a narrow HOMO-LUMO gap and a substantial DRC. In most cases, the low ionization potential is related to a large DRC, but the ionization potential alone cannot explain the variations in the DRC. The electronegativity essentially follows the profile of the ionization potential, while the electrophilicity matches softness.
For clarification of the laboratory stability of quinone methides, analysis of the local reactivity indices was needed in addition to the DRC. Taking into account the knowledge about the characteristic reaction mechanisms for this class of compounds, most important should be the Fukui indices for nucleophilic (water) and radical (dimerization) attacks. As anticipated, in the compounds without substituents in the exocyclic methylene group (the a series), the carbon atom in the latter group is the most reactive. In general, the a series features higher values of the local Fukui indices than the substituted analogs (b series). The comparison of DRC with f + max outlines the following extreme cases: (1) low DRC (0.114) + high f + max (0.160) = instability (Ia) (2) low DRC (0.084) + medium-to-high f + max (0.120) = moderate stability (XII) (3) low DRC (0.074) + low f + max (0.046) = stability (XIb) (4) medium DRC (0.394) + medium-to-high f + max (0.125) = instability (VIa) (5) medium DRC (0.508) + low f + max (0.046) = moderate stability (VIb).
As the Fukui indices for the three types of attack display analogous behavior within each series, the same relations apply for all three. The experimentally observed laboratory stability is represented by the color code. Highly stable quinone methides have P > 55 and they were obtained experimentally in the crystal phase. Compounds with 30 < P < 55 have moderate stability, and those with P < 30 are unstable.
For clarification of the laboratory stability of quinone methides, analysis of the local reactivity indices was needed in addition to the DRC. Taking into account the knowledge about the characteristic reaction mechanisms for this class of compounds, most important should be the Fukui indices for nucleophilic (water) and radical (dimerization) attacks. As anticipated, in the compounds without substituents in the exocyclic methylene group (the a series), the carbon atom in the latter group is the most reactive. In general, the a series features higher values of the local Fukui indices than the substituted analogs (b series). The comparison of DRC with f + max outlines the following extreme cases: (1) low DRC (0.114) + high f + max (0.160) = instability (Ia) (2) low DRC (0.084) + medium-to-high f + max (0.120) = moderate stability (XII) (3) low DRC (0.074) + low f + max (0.046) = stability (XIb) (4) medium DRC (0.394) + medium-to-high f + max (0.125) = instability (VIa) (5) medium DRC (0.508) + low f + max (0.046) = moderate stability (VIb).
As the Fukui indices for the three types of attack display analogous behavior within each series, the same relations apply for all three.
The analysis reveals a fine balance between the DRC and reactivity, which may result in both stable and unstable compounds with an open-shell character. To further examine the complex interrelations within each series, we performed a chemometrics itemization of the data employing cluster analysis. A table with all descriptors and compounds subject to the cluster analysis is presented in the SM, Table S1.
The chemometrics analysis results are presented in Figure 4 (as well as in the SM, Figures S1-S2, Tables S2-S4) and show a very good separation of the quinone methides based on their electronic properties, reactivity, and stability. Three groups emerge: Cluster 1, Cluster 2, and Cluster 3. Cluster 3 encompasses compounds Ia, IIa, IIIa, IVa, VIIa, Xia, and XII. This is the group of the low-to-moderate stable molecules. All compounds in this cluster belong to the a series and contain a highly reactive non-substituted exocyclic carbon. They feature lowest electron affinity, highest ionization potential (prominent reducers), and, hence, lowest softness and electrophilicity. Their averaged diradical character amounts to 0.11 and is the lowest for the series range. However, they are also characterized by the largest atomic Fukui indices and are therefore highly reactive. The compounds from Cluster 3 are not interesting for application because they do not possess a significant diradical character or laboratory stability. compounds subject to the cluster analysis is presented in the SM, Table S1.
The chemometrics analysis results are presented in Figure 4 (as well as in the SM, Figures S1-S2, Tables S2-S4) and show a very good separation of the quinone methides based on their electronic properties, reactivity, and stability. Three groups emerge: Cluster 1, Cluster 2, and Cluster 3. Cluster 3 encompasses compounds Ia, IIa, IIIa, IVa, VIIa, Xia, and XII. This is the group of the low-to-moderate stable molecules. All compounds in this cluster belong to the a series and contain a highly reactive non-substituted exocyclic carbon. They feature lowest electron affinity, highest ionization potential (prominent reducers), and, hence, lowest softness and electrophilicity. Their averaged diradical character amounts to 0.11 and is the lowest for the series range. However, they are also characterized by the largest atomic Fukui indices and are therefore highly reactive. The compounds from Cluster 3 are not interesting for application because they do not possess a significant diradical character or laboratory stability. Cluster 2 contains the moderate-to-high stability compounds and comprises all molecules from the b series: Ib, IIb, IIIb, IVb, VIIb, and XIb. All these quinone methides have structural analogues from Cluster 3. The average diradical character of the compounds in Cluster 2 is 0.18 (i.e., higher than that for Cluster 3), but their lower reactivity, and hence higher stability, are ensured by phenyl functionalization of the methylene carbon. The latter is obvious from the lowest values of the atomic Fukui indices, too. The molecules from Cluster 2 are also characterized by the lowest ionization potential and moderate electron affinity, which results in moderate-to-high softness in the series. Clus- Cluster 2 contains the moderate-to-high stability compounds and comprises all molecules from the b series: Ib, IIb, IIIb, IVb, VIIb, and XIb. All these quinone methides have structural analogues from Cluster 3. The average diradical character of the compounds in Cluster 2 is 0.18 (i.e., higher than that for Cluster 3), but their lower reactivity, and hence higher stability, are ensured by phenyl functionalization of the methylene carbon. The latter is obvious from the lowest values of the atomic Fukui indices, too. The molecules from Cluster 2 are also characterized by the lowest ionization potential and moderate electron affinity, which results in moderate-to-high softness in the series. Cluster 2 is a good example of a reasonable compromise between a relatively significant diradical character and stability. This result clearly demonstrates that stable molecular materials with DRC exist, and that the combination between quantum-chemical and chemometrics approaches is a powerful tactic to find them. Finally, compounds Va, Vb, VIb, VIa, VIb, VIII, IX, X, and XIII belong to Cluster 1, which differentiates the group of the moderately-stabileto-unstable molecules. In the whole series, these compounds have the highest electron affinities (prominent oxidizers), lowest hardness, highest softness and electrophilicity. Their averaged DRC is intermediate, 0.38, which is also a high value for the series. Therefore, although attractive from a properties perspective, such compounds are useless in practice.
From the discussion on the quantum-chemical results, it is clear that the used descriptors are correlated. This correlation, however, does not have a negative effect on the clustering (hierarchical or K-means) since cluster labels are assigned on the basis of the Euclidean distance between the points in the clusters and their corresponding centroids in the 15-dimensional feature space.
Inspection of the ground state resonance structures discloses that the distinct partition of quinone methides in terms of electronic properties and laboratory stability can be interpreted by implementing Clar's rule [34,41] for the aromatic sextets ( Figure 5). the highest electron affinities (prominent oxidizers), lowest hardness, highest softness and electrophilicity. Their averaged DRC is intermediate, 0.38, which is also a high value for the series. Therefore, although attractive from a properties perspective, such compounds are useless in practice.
From the discussion on the quantum-chemical results, it is clear that the used descriptors are correlated. This correlation, however, does not have a negative effect on the clustering (hierarchical or K-means) since cluster labels are assigned on the basis of the Euclidean distance between the points in the clusters and their corresponding centroids in the 15-dimensional feature space.
Inspection of the ground state resonance structures discloses that the distinct partition of quinone methides in terms of electronic properties and laboratory stability can be interpreted by implementing Clar's rule [34,41] for the aromatic sextets ( Figure 5). Figure 5. Kekulé structures of quinone methides containing Clar sextets (shaded) having a stabilizing effect on the ground state. The experimentally observed laboratory stability of the ground state is represented a color code and the Clar sextets are shaded correspondingly. All bi-and tricyclic compounds without Clar sextets in the basic structure belong to Cluster 1, the group of the unstable molecules. The rest have at least one Clar sextet in the quinone methide fragment and fall in Cluster 2 and 3-the groups of stable and moderately stable molecules.
All stable and moderately stable compounds (Cluster 2 and 3) possess at least one Clar sextet in the main quinone methide structure (excluding the phenyl substituents). Ib and IIb are an exception, but in these, the stabilization originates from the two phenyl substituents that import two Clar sextets, thus greatly stabilizing the single six-membered quinoid ring of the main structure. For the remaining compounds, the introduction of two phenyl substituents results in no more than moderate stability. Our Figure 5. Kekulé structures of quinone methides containing Clar sextets (shaded) having a stabilizing effect on the ground state. The experimentally observed laboratory stability of the ground state is represented a color code and the Clar sextets are shaded correspondingly. All bi-and tricyclic compounds without Clar sextets in the basic structure belong to Cluster 1, the group of the unstable molecules. The rest have at least one Clar sextet in the quinone methide fragment and fall in Cluster 2 and 3-the groups of stable and moderately stable molecules.
All stable and moderately stable compounds (Cluster 2 and 3) possess at least one Clar sextet in the main quinone methide structure (excluding the phenyl substituents). Ib and IIb are an exception, but in these, the stabilization originates from the two phenyl substituents that import two Clar sextets, thus greatly stabilizing the single six-membered quinoid ring of the main structure. For the remaining compounds, the introduction of two phenyl substituents results in no more than moderate stability. Our recent study on singlet fission chromophores also demonstrates the importance of Clar's sextet rule for the design of photovoltaics with diradical character [16]. Therefore, it is worth investing in the development of high-throughput approaches for counting Clar's sextets and their implementation in qualitative structure-property relationships models for finding new solar cell materials. To the best of our knowledge, currently, computational approaches for counting Clar's sextets exist but are limited to polycyclic aromatic hydrocarbons [42].

Conclusions
The laboratory stability of quinone methides was rationalized using their diradical character, global and local reactivity indices, chemometrics analysis, and Clar's rule. All investigated compounds possess non-zero DRCs, and the stable or moderately stable among them possess at least one Clar sextet in the resonance structure. It was established that the average DRC value of 0.18 is in accordance with laboratory stability, which demonstrates that a compromise between a non-zero diradical character and real application can be found. For this purpose, the theoretical investigations must include global and local reactivity indices as this class of descriptors is related to the DRC. The study reveals that the reactivity indices and Clar's sextet rule open new horizons for high-throughput materials design. They can be implemented in chemometrics approaches for the quest of stable materials for applications in singlet fission photovoltaics and photonics.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/sym13081448/s1, Table S1: Dataset used for the K-means clustering of quinone methides; Tables S2-S4: Members distribution in the obtained clusters from K-means clustering; Figure S1: Hierarchical clustering for the variable-dendrogram representation; Figure S2: Hierarchical clustering of the objects-a dendrogram representation.