Machine-learning-assisted performance improvements for multi-resonance thermally activated delayed fluorescence molecules

Abstract

With favorable colour purity, multi-resonance thermally activated delayed fluorescence (MR-TADF) molecules exhibit enormous potential in high-definition displays. Due to the relatively small chemical space of MR-TADF molecules, it is challenging to improve molecular performance through domain-specific expertise alone. To address this problem, we focused on optimizing the classic molecule, DABNA-1, using machine learning (ML). Molecular morphing operations were initially employed to generate the adjacent chemical space of DABNA-1. Subsequently, a machine learning model was trained with a limited database and used to predict the properties throughout the generated chemical space. It was confirmed that the top 100 molecules suggested by machine learning present excellent electronic structures, characterized by small reorganization energy and singlet–triplet energy gaps. Our results indicate that the improvement in electronic structures can be elucidated through the view of the molecular orbital (MO). The results also reveal that the top 5 molecules present weaker vibronic peaks of the emission spectrum, demonstrating higher colour purity when compared to DABNA-1. Notably, the M2 molecule presents a high RISC rate, indicating its promising future as a high-efficiency MR-TADF molecule. Our machine-learning-assisted approach facilitates the rapid optimization of classical molecules, addressing a crucial requirement within the organic optoelectronic materials community.

---

1. Introduction

Multi-resonance thermally activated delayed fluorescence (MR-TADF) molecules hold great promise in addressing the broad emission issue commonly encountered by conventional TADF molecules.^1–4 This is closely associated with the weak vibronic coupling that arises from local electron density reorganization, a key feature of the MR-TADF mechanism.³ As a result, MR-TADF molecules are emerging as an attractive choice for high-definition organic light-emitting diode (OLED) emitters. In general, an excellent MR-TADF molecule should present weak vibronic coupling between the lowest singlet excited state (S₁) and the ground state (S₀), as well as a small energy gap (ΔE_ST) between S₁ and the lowest triplet excited state (T₁).^2,5 The weak vibronic coupling is associated with narrowband emission, thereby ensuring high colour purity.^6–8 A small ΔE_ST can accelerate the reverse intersystem crossing (RISC) process of triplet excitons upconverting to singlet states, which reduces triplet exciton concentrations and mitigates efficiency roll-off.^3,9–11

In comparison to the conventional TADF molecules with donor–acceptor structures, the chemical space explored for MR-TADF molecules (the number of known MR-TADF structures) remains limited.^3,12 The main design strategy for MR-TADF molecules is to combine heteroatoms with complementary electronegativity, currently with the commonly available combinations of B/N, B/O, and C=O/N.^1,3 As a result, it will be challenging to further improve the colour purity and the RISC rate through expert experience. As a classic well-performing MR-TADF molecule, B/N-doped compound DABNA-1 (Fig. 1) has experimentally exhibited a narrow full width at half maximum (FWHM) of 190 meV and a remarkable RISC rate of 9.9 × 10³ s⁻¹.¹³

Fig. 1 Schematic overview of machine-learning-assisted performance improvements for MR-TADF molecules.

Quantum chemical calculations are perceived as a more cost-efficient alternative to experimental methods. However, large-scale computational screening costs remain substantial, particularly regarding the vibronic coupling and ΔE_ST calculations for MR-TADF molecules. Here, the strength of vibronic coupling can be characterized by the reorganization energy.^6,7 The measure of the vibronic coupling strength (the reorganization energy) requires the calculation of excited state forces.⁶ The evaluation of ΔE_ST should be conducted by using the wave function-based methods due to the double-excitation character of electronic transitions.^4,14 As a promising approach to address the challenge posed by high-cost calculations and limited computational resources, machine learning (ML) has demonstrated tremendous potential as a valuable tool in the fields of medicinal chemistry, materials chemistry, and organic synthesis.^10,15–17 Over the past few years, it has received significant attention and exploration.

In this work, a cost-effective approach is presented to improve the performance of MR-TADF molecules using ML, with DABNA-1 as a representative example. The schematic overview is shown in Fig. 1. First, an ML model was trained with a limited database and subsequently used to predict the properties of the entire adjacent chemical space of DABNA-1. Second, a series of structure–property relationship analyses were conducted for the screened excellent molecules from the perspective of the molecular orbital (MO). Finally, the emission spectra and ISC/RISC rates of DABNA-1 and the top 5 MR-TADF molecules were calculated to validate the effectiveness and practicality of our machine-learning-assisted performance improvement approach.

2. Results and discussion

2.1 Exploring the adjacent chemical space of DABNA-1 using machine learning with molecular graph neural networks

An MR-TADF molecule can generally be fragmented into a molecular backbone (core) and side groups (attached to the core), where the core usually dominates its photophysical properties. To efficiently generate the adjacent chemical space of DABNA-1, we iteratively applied a set of molecular construction rules (Fig. 2), i.e., morphing operations, to its core Core_0 (Fig. 1 and 2). These morphing operations were considered to use the resonance effect of heteroatoms to separate the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO) in the modification of molecular structures. These operations were applied to all possible atomic sites or fragments of the molecules. The generation process was initiated with Core_0 and carried out with 10 cycles of the morphing operations in Fig. 2. This process produced up to 301 852 molecules within the adjacent chemical space of DABNA-1. It should be noted that all these molecules were constrained with molecular symmetry, as symmetrical molecules always present weak vibronic coupling and are conducive to synthesis.^17–19

Fig. 2 Molecular morphing operations applied for generating the adjacent chemical space of DABNA-1. Here, the blue colour denotes the molecule before and after structure modifications, and the red colour denotes molecular morphing operations.

Following the establishment of the adjacent chemical space of DABNA-1, 5000 molecules were randomly selected to construct a database for machine learning. It is commonly believed that molecules with strong structural rigidity usually present a remarkable similarity between their S₀ and S₁ potential energy surfaces. As a result, the reorganization energy for emission (λ_em) is nearly equal to the reorganization energy for absorption (λ_abs) (Fig. 3a).^6,20,21 To eliminate the high-cost computational expenses of geometries and vibrational frequencies of the excited state, the reorganization energy under the vertical gradient (VG) approximation at S₀ geometry (λ^VG) as well as the vertical S₁–T₁ energy gap (ΔE^vert_ST) at S₀ geometry (Fig. 3a) were used as descriptors to characterize colour purity and RISC rates, respectively.⁶ Additionally, the two individual properties (λ^VG and ΔE^vert_ST) were converted into two single scores ranging from zero to one by using a logistic function (Fig. 3b and d).²² A higher score denotes a better individual property, with a score of above 0.8 indicating exceptional performance. Unfortunately, the calculated score of the two investigated properties of most molecules is below 0.8 (Fig. 3c and e). Here, to reflect the overall excellence of the two parameters of a molecule, the multi-parameter optimization (MPO) score was calculated using the formula , where a_i is the score for the individual property and n is the number of individual properties.²² A high MPO score indicates that all the individual properties of the molecule are favorable.

Fig. 3 (a) The potential energy diagram of the S₀, S₁, and T₁ states. The desirability score for (b) λ^VG and (d) ΔE^vert_ST. The score distribution for (c) λ^VG and (e) ΔE^vert_ST in the database with 5000 molecules.

Furthermore, a molecular graph convolutional network (GCN) was trained by using the database with 5000 molecules, incorporating relevant information such as the atomic number and hybridization, to generate a prediction model.²³ The coefficient of determination and the mean absolute error were used to evaluate the GCN model, where E_ML and E_QC denote the model-predicted value and actual value in the test data, respectively, and N represents the number of molecules in the test data.

As shown in Fig. 4a and b, the R² values of λ^VG and ΔE^vert_ST are calculated to be 0.77 and 0.72, respectively, while the MAEs of λ^VG and ΔE^vert_ST are 18 meV and 50 meV, respectively. These results indicate that the prediction model can effectively estimate the two molecular properties. Applying this prediction model throughout the adjacent chemical space of DABNA-1, the top 100 molecules with the highest MPO scores were selected to conduct further quantum chemical validation. The data distributions of the initial 5000 molecules and the predicted top 100 molecules are depicted in Fig. 4c. The MPO scores of the predicted top 100 molecules are around 0.8, which demonstrates that our machine-learning-assisted performance improvement approach is very promising.

Fig. 4 Comparison of (a) λ^VG and (b) ΔE^vert_ST obtained from quantum chemistry calculations and GCN models. (c) MPO score distributions of the database with 5000 molecules and the predicted top 100 molecules. The data of the top 100 molecules predicted by the GCN model were recalculated using quantum chemistry methods.

2.2 Structure–property analysis of the top 100 MR-TADF molecules

To clarify the structural characteristics of these high-performance MR-TADF molecules, the Morgan circular fingerprints of the top 100 molecules were calculated to compare the structural similarities (Fig. 5a).²⁴ The heatmap analysis shows that these molecular structures can be grouped into three distinct clusters due to the high within-group similarity. This indicates that these molecules are mainly based on three cores, with Core_1 having been identified as a particularly effective MR-TADF core structure.²⁵ To explore the structure–property relationship of these cores, the geometries of singlet and triplet excited states were further optimized to calculate the vibronic coupling strength for emission (λ_em) and the adiabatic S₁–T₁ energy gap (ΔE^ad_ST).

Fig. 5 (a) Similarity heatmap for the top 100 MR-TADF molecules. (b) The main core structures for the top 100 molecules.

As shown in Fig. 6b, Core_1 and Core_3 present a noticeable decrease in reorganization energy (λ_em) compared to Core_0. This decrease mainly comes from the bond stretching modes around 1400–1700 cm⁻¹ (Fig. S4, ESI†). To deeply understand the reorganization energy, λ_em was further decomposed into internal coordinates.^26,27 In the internal coordinate representation, λ_em can be described as the result of alterations in bond lengths, bond angles, and dihedral angles.²⁷ The reorganization energy resulting from alterations in bond lengths (λ_BL) significantly contributes to the decreased λ_em from Core_0 to Core_1 and Core_3 (Fig. 6c). Particularly, this decrease in λ_BL can be primarily attributed to the alterations in the C6–C13 and C10–C11 bonds (Fig. 6d, g, and m). It is worth noting that the bond length is associated with bond order. The bond order between the μth and νth atoms can be written as , where n_i represents the occupation number of the ith MO, Cⁱ_μ and Cⁱ_ν denote the MO coefficients.^20,27 Therefore, λ_BL can be rationalized from the perspective of MO theory, specifically by using λ_BL ∝ (C^L_μC^L_ν − C^H_μC^H_ν)², where H and L denote the HOMO and LUMO, respectively.^20,27 As shown in Fig. 6a, the C6–C13 and C10–C11 bonds show the strong bonding character in the HOMO, while showing the weak bonding character in the LUMO. As molecular conjugation increases from Core_0 to Core_1 and Core_3, there is a noticeable decrease in the HOMO distribution located on the C11 and C13 atoms (Fig. 6e and k). As a result, the C6–C13 and C10–C11 bonds present a decreased bond order difference (associated with bond length alteration) between the S₀ and S₁ states, thus showing a decrease in λ_BL.

Fig. 6 (a) HOMO and LUMO distributions (iso_value = 0.03 au) of the four cores. (b) Calculated λ_em and ΔE^ad_ST. (c) The contribution of the bond length, bond angle and dihedral angle to λ_em. (d), (g), (j), and (m) The contribution of each chemical bonds to λ_BL. (e), (f), (h), (i), (k), and (l) The difference of the orbital composition (ΔOC_{Core_1/Core_2/Core_3} – _{Core_0}) on atoms.

The ΔE^ad_ST values of Core_0, Core_1, Core_2 and Core_3 are calculated to be 198, 124, 214, and 141 meV, respectively. Given that the spatial overlap between the HOMO and LUMO determines the energy gap between the S₁ and T₁ states.²⁸ The squares of the overlap integral of the HOMO and the LUMO (OverlapHL2) for Core_0, Core_1, Core_2 and Core_3 are calculated to be 13.40 × 10⁻⁴, 7.11 × 10⁻⁴, 9.03 × 10⁻⁴, and 7.07 × 10⁻⁴, respectively. The reduced HOMO distributions on the N1, N9, C11, C13, C17, and C18 atoms, as well as the reduced LUMO distributions on the C17 and C18 atoms (Fig. 6e, f, k, and l), are the main contributors to the decrease in ΔE^ad_ST from Core_0 to Core_1 and Core_3.

It should be noted that it is not easy to simultaneously weaken vibronic coupling and decrease ΔE_ST for the conventional TADF molecules. This is because decreasing ΔE_ST in conventional TADF molecules is achieved by enhancing the long-distance spatial separation of the HOMO and LUMO. However, this long-range electron density reorganization is accompanied by the strong vibronic coupling (the large reorganization energy).^2,3

2.3 Photophysical properties of the top 5 MR-TADF molecules

The top 5 MR-TADF molecules are compared to the three cores to examine their improvement in two significant physical properties. Specifically, M1, M2, and M3 are based on Core_1, M4 is based on Core_2, and M5 is based on Core_3. As the degree of conjugation increases, M1, M2, M3, and M4 present smaller λ_em and lower ΔE^ad_ST compared to their respective corresponding parent structures, while M5 has no evident decrease in λ_em and ΔE^ad_ST compared to Core_3 (Fig. 6b and 7b). Here, the Overlap_HL² for M1, M2, M3, M4 and M5 are calculated to be 5.24 × 10⁻⁴, 4.79 × 10⁻⁴, 5.37 × 10⁻⁴, 7.56 × 10⁻⁴, and 6.02 × 10⁻⁴, respectively. Consistent with the analysis in Section 2.2, the reduced HOMO distribution located on the C11 and C13 atoms (Fig. S7, ESI†) serves as a critical determinant for the differences in λ_em and ΔE^ad_ST between M1, M2, M3, and M4 and their respective corresponding parent structures (Fig. 6b and 7b). Because this reduced HOMO distribution weakens the difference in bond order (associated with λ_BL) of the C6–C13 and C10–C11 bonds (Fig. S6, ESI†) and decreases the HOMO–LUMO overlap (associated with ΔE^ad_ST).

Fig. 7 (a) Chemical structures of the top 5 MR-TADF molecules. (b) Calculated λ_em and ΔE^ad_ST. (c) Calculated fluorescence emission spectra.

To visually verify the superior performance of MR-TADF molecules suggested by our machine-learning-assisted approach, the fluorescence emission spectra and ISC/RISC processes from T₁ to S₁ were studied. With the smaller λ_em, the top 5 MR-TADF molecules show the weaker vibronic peak (the narrower FWHM) compared to DABNA-1 (Fig. 7c and Table 1), especially for the bottom of the low energy region in emission spectra. As a result, the top 5 MR-TADF molecules show high colour purity. Meantime, based on the smaller ΔE^ad_ST, the 5 MR-TADF molecules present an increased ratio of RISC rate k_RISC(T₁ → S₁) to ISC rate k_ISC(S₁ → T₁) compared to DABNA-1 (Table 1), which is beneficial to enhance the exciton utilization. In particular, M2 demonstrates a promising future as a high-efficiency MR-TADF molecule. Owing to its ΔE^ad_ST of 50 meV, the calculated k_RISC(T₁ → S₁) value of M2 can reach 8.97 × 10⁵ s⁻¹, indicating great potential to avoid efficiency roll-off. It should be noted that the top 5 molecules are demonstrated with high radiative rates (Table 1), although our primary focus during the ML process was not specifically on achieving high radiative rates.

Table 1 Calculated FWHM, ISC rate k_ISC(S₁ → T₁), RISC rate k_RISC(T₁ → S₁), and radiative rate k_r(S₁ → S₀)

	DABNA-1	M1	M2	M3	M4	M5
FWHM (meV)	196	161	148	130	134	102
k ISC(S1 → T1) (s^−1)	3.52 × 10^5	5.39 × 10^4	3.83 × 10^6	9.34 × 10^3	1.71 × 10^5	2.56 × 10^4
k RISC(T1 → S1) (s^−1)	6.64 × 10^2	3.65 × 10^3	8.97 × 10^5	3.25 × 10^3	1.76 × 10^3	3.30 × 10^2
k r(S1 → S0) (s^−1)	6.72 × 10^7	1.53 × 10^8	1.57 × 10^8	1.77 × 10^8	1.07 × 10^8	4.37 × 10^8

---

3. Conclusions

To summarize, machine learning was employed to advance the performance of MR-TADF molecules, with a specific focus on improving the colour purity and RISC rate simultaneously. To quantify the two critical properties, λ^VG and ΔE^vert_ST were considered as descriptors to characterize the colour purity and RISC rate, respectively. Using the molecular morphing operations, the adjacent chemical space of DABNA-1 was generated up to 301 852 molecules. At an acceptable cost, 5000 molecules were randomly selected from the adjacent space to perform molecular graph neural network learning and generate a prediction model. The prediction model was then used to estimate the molecular properties throughout the adjacent chemical space, achieving R² of 0.77 and 0.72 for λ^VG and ΔE^vert_ST, respectively, with the corresponding MAE of 18 and 50 meV. Moreover, the MPO scores of the top 100 molecules suggested by our GCN model are around 0.8, indicating the effectiveness of our approach in identifying high-efficiency MR-TADF molecules.

Our structure–property analysis revealed that the top 100 MR-TADF molecules are mainly based on three core structures, and it is possible to elucidate the reorganization energy and S₁–T₁ energy gaps through the view of MO theory. The reduced HOMO distribution located on the C11 and C13 atoms serves as a critical determinant for the decrease in λ_em and ΔE^ad_ST from Core_0 to Core_1 and Core_3. Moreover, this reduced HOMO distribution is associated with the decrease in λ_em and ΔE^ad_ST from Core_1 to M1, M2, and M3, as well as from Core_2 to M4. On the other hand, the top 5 MR-TADF molecules show a weaker vibronic peak in the emission spectrum compared to DABNA-1. Notably, M2 with a calculated k_RISC(T₁ → S₁) of 8.97 × 10⁵ s⁻¹ offers an effective solution to efficiency roll-off, thus signifying its potential as a high-efficiency MR-TADF molecule. By integrating molecular structure modification, ML, and insights from MO, we can improve the performance of MR-TADF molecules and gain a deeper understanding of the underlying improvement process.

In the era of artificial intelligence, machine learning plays an important role in discovering novel materials. Numerous efforts have been made to improve the model performance and develop low-cost and effective first-principles descriptors. However, to advance the performance of specific materials, it is crucial to generate the chemical space with high-efficiency molecules through specific morphing operations. It is worth mentioning that an efficient structure–property relationship analysis can guide the development of specific morphing operations.

4. Methods

4.1 Molecular structure and electronic structure calculations

The morphing operations used to generate adjacent chemical space were conducted using the “Reaction SMARTS” method in the RDKit program.^17,29 Molecules in this space were represented by the simplified molecular input line entry system (SMILES). With the SMILES string selected to construct the database, the initial molecular geometry was generated by using the MMFF94 force field, and then optimized at the PBE0-D3BJ/def2-SVP level in the Gaussian 16 program.^30–35 The reorganization energy under the vertical gradient approximation was evaluated as , where the displacement of the kth mode was calculated by using the gradient of the S₁ state at S₀ geometry (Q₀), .^20,36 The ΔE^vert_ST at S₀ geometry was calculated by using the spin-component scaled second-order approximate coupled cluster (SCS-CC2) method with the cc-pVDZ basis set in the MRCC program.^37–39

To perform detailed photophysical property analyses, the S₁ and T₁ states of DABNA-1, the four cores, and the top 5 MR-TADF molecules were respectively optimized by using the time-dependent density functional theory (TDDFT) and unrestricted DFT in the Gaussian 16 program. The S₁ and T₁ excitation energies were evaluated at the SCS-CC2/cc-pVTZ level in the MRCC program.^38,39 The spin–orbit coupling (SOC) matrix elements between the S₁ and T₁ states were calculated at the PBE0-D3BJ/def2-SVP level by using the PySOC program.⁴⁰ The Overlap_HL² was obtained from the Multiwfn program.⁴¹

4.2 Photophysical property calculations

Based on Fermi's golden rule, the emission spectrum can be written as⁴²where μ_fi and P_iν denote the transition dipole moment and the Boltzmann distribution, respectively. |〈Θ_fν′|Θ_iν〉|² is the Franck–Condon factor between the S₁ and S₀ states, which reflects the strength of vibronic coupling and determines the intensity of the vibronic peak. Based on the electronic structure calculation data, the emission spectrum calculation was carried out using the MOMAP program.⁴³

The radiative decay rate (k_r(S₁ → S₀)) can be obtained from the Einstein spontaneous emission equation⁴⁴

where f_osc is the oscillator strength and ΔE is the vertical excitation energy in cm⁻¹ at the S₁ geometry.

The ISC/RISC rate can be calculated as follows:^45,46

where H^SO_fi denotes the spin–orbit coupling matrix element. ρ_fi (t,T) is the Franck–Condon factor between the S₁ and T₁ states.

4.3 Molecular graph neural network

The graph convolutional network training process was performed by using the Deep Graph Library (DGL) package based on PyTorch.^23,47 The propagation process of this network is written as , which means a linear projection (W^l) of the aggregated information followed by a non-linear activation function (σ). Ã = A + I is the sum of the adjacency matrix A and the identity matrix I. D̃ is the diagonal node degree matrix of Ã. H^l and W^l are the lth neural network layer and its learnable parameters, respectively. In the case of the input layer, H^l corresponds to the molecular node information. Here, the hyperparameter search, including the determination of the number and size of hidden layers, was performed by using Bayesian optimization.⁴⁸ More details and code are available at https://github.com/Wanlin-Cai/ML_GCN.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work is supported by the National Natural Science Foundation of China (22032004, 22021001, and 21533006).

References

1. S. Madayanad Suresh, D. Hall, D. Beljonne, Y. Olivier and E. Zysman-Colman, Adv. Funct. Mater., 2020, 30, 1908677.
2. J. M. Ha, S. H. Hur, A. Pathak, J.-E. Jeong and H. Y. Woo, NPG Asia Mater., 2021, 13, 1–36.
3. D. Hall, K. Stavrou, E. Duda, A. Danos, S. Bagnich, S. Warriner, A. M. Slawin, D. Beljonne, A. Köhler and A. Monkman, Mater. Horiz., 2022, 9, 1068–1080.
4. A. Pershin, D. Hall, V. Lemaur, J.-C. Sancho-Garcia, L. Muccioli, E. Zysman-Colman, D. Beljonne and Y. Olivier, Nat. Commun., 2019, 10, 597.
5. V. V. Patil, H. L. Lee, I. Kim, K. H. Lee, W. J. Chung, J. Kim, S. Park, H. Choi, W. J. Son and S. O. Jeon, Adv. Sci., 2021, 8, 2101137.
6. S. A. Ahmad, J. Eng and T. J. Penfold, J. Mater. Chem. C, 2022, 10, 4785–4794.
7. W. L. Cai, A. C. Zhao, K. Ren, R. X. He, M. Li and W. Shen, J. Phys. Chem. C, 2019, 123, 17968–17975.
8. Z. Pei, Q. Ou, Y. Z. Mao, J. J. Yang, A. L. D. L. Lande, F. Plasser, W. Z. Liang, Z. G. Shuai and Y. H. Shao, J. Phys. Chem. Lett., 2021, 12, 2712–2720.
9. P. K. Samanta, D. Kim, V. Coropceanu and J.-L. Brédas, J. Am. Chem. Soc., 2017, 139, 4042–4051.
10. R. Gómez-Bombarelli, J. Aguilera-Iparraguirre, T. D. Hirzel, D. Duvenaud, D. Maclaurin, M. A. Blood-Forsythe, H. S. Chae, M. Einzinger, D.-G. Ha and T. Wu, et al. , Nat. Mater., 2016, 15, 1120–1127.
11. R. Pollice, P. Friederich, C. Lavigne, G. dos Passos Gomes and A. Aspuru-Guzik, Matter, 2021, 4, 1654–1682.
12. X. Zeng, X. Wang, Y. W. Zhang, G. Y. Meng, J. B. Wei, Z. Y. Liu, X. Q. Jia, G. M. Li, L. Duan and D. D. Zhang, Angew. Chem., Int. Ed., 2022, 134, e202117181.
13. T. Hatakeyama, K. Shiren, K. Nakajima, S. Nomura, S. Nakatsuka, K. Kinoshita, J. Ni, Y. Ono and T. Ikuta, Adv. Mater., 2016, 28, 2777–2781.
14. S. M. Pratik, V. Coropceanu and J.-L. Brédas, ACS Mater. Lett., 2022, 4, 440–447.
15. S. D. Xu, J. L. Li, P. F. Cai, X. L. Liu, B. Liu and X. N. Wang, J. Am. Chem. Soc., 2021, 143, 19769–19777.
16. M. Sumita, K. Terayama, N. Suzuki, S. Ishihara, R. Tamura, M. K. Chahal, D. T. Payne, K. Yoshizoe and K. Tsuda, Sci. Adv., 2022, 8, eabj3906.
17. C. Kunkel, J. T. Margraf, K. Chen, H. Oberhofer and K. Reuter, Nat. Commun., 2021, 12, 2422.
18. C.-C. Wu, E. Y. Li and P.-T. Chou, Chem. Sci., 2022, 13, 7181–7189.
19. M. Uejima, T. Sato, K. Tanaka and H. Kaji, Phys. Chem. Chem. Phys., 2013, 15, 14006–14016.
20. W.-C. Chen and Y.-C. Cheng, J. Phys. Chem. A, 2020, 124, 7644–7657.
21. C. Adamo and D. Jacquemin, Chem. Soc. Rev., 2013, 42, 845–856.
22. H. S. Kwak, Y. An, D. J. Giesen, T. F. Hughes, C. T. Brown, K. Leswing, H. Abroshan and M. D. Halls, Front. Chem., 2022, 9, 800370.
23. M. Li, J. Zhou, J. Hu, W. Fan, Y. Zhang, Y. Gu and G. Karypis, ACS Omega, 2021, 6, 27233–27238.
24. H. L. Morgan, J. Chem. Doc., 1965, 5, 107–113.
25. Y. W. Zhang, D. D. Zhang, J. B. Wei, Z. Y. Liu, Y. Lu and L. Duan, Angew. Chem., Int. Ed., 2019, 58, 16912–16917.
26. H. Geng, Y. L. Niu, Q. Peng, Z. G. Shuai, V. Coropceanu and J.-L. Brédas, J. Chem. Phys., 2011, 135, 104703.
27. W. L. Cai, C. Zhong and D. Y. Wu, Mater. Chem. Front., 2023, 7, 3762–3773.
28. Z. Y. Yang, Z. Mao, Z. L. Xie, Y. Zhang, S. W. Liu, J. Zhao, J. R. Xu, Z. G. Chi and M. P. Aldred, Chem. Soc. Rev., 2017, 46, 915–1016.
29. G. Landrum, RDKit: Open-Source Cheminformatics Software (accessed September, 2021) https://www.rdkit.org.
30. T. A. Halgren, J. Comput. Chem., 1996, 17, 490–519.
31. C. Adamo and V. Barone, J. Chem. Phys., 1999, 110, 6158–6170.
32. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
33. S. Grimme, S. Ehrlich and L. Goerigk, J. Comput. Chem., 2011, 32, 1456–1465.
34. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
35. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, et al., Gaussian 16 Rev. A.03, Gaussian Inc., Wallingford CT, 2016.
36. F. Santoro and D. Jacquemin, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2016, 6, 460–486.
37. M. Kállay, P. R. Nagy, D. Mester, Z. Rolik, G. Samu, J. Csontos, J. Csóka, P. B. Szabó, L. Gyevi-Nagy and B. Hégely, J. Chem. Phys., 2020, 152, 074107.
38. O. Christiansen, H. Koch and P. Jørgensen, Chem. Phys. Lett., 1995, 243, 409–418.
39. T. H. Dunning Jr, J. Chem. Phys., 1989, 90, 1007–1023.
40. X. Gao, S. M. Bai, D. Fazzi, T. Niehaus, M. Barbatti and W. Thiel, J. Chem. Theory Comput., 2017, 13, 515–524.
41. T. Lu and F. Chen, J. Comput. Chem., 2012, 33, 580–592.
42. S. H. Lin, C. H. Chang, K. K. Liang, R. Chang, Y. J. Shiu, J. M. Zhang, T. S. Yang, M. Hayashi and F. C. Hsu, Adv. Chem. Phys., 2002, 1–88.
43. Y. L. Niu, W. Q. Li, Q. Peng, H. Geng, Y. P. Yi, L. J. Wang, G. J. Nan, D. Wang and Z. G. Shuai, Mol. Phys., 2018, 116, 1078–1090.
44. R. C. Hilborn, Am. J. Phys., 1982, 50, 982–986.
45. Q. Peng, Q. H. Shi, Y. L. Niu, Y. P. Yi, S. R. Sun, W. Q. Li and Z. G. Shuai, J. Mater. Chem. C, 2016, 4, 6829–6838.
46. Z. G. Shuai and Q. Peng, Natl. Sci. Rev., 2017, 4, 224–239.
47. A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein and L. Antiga, Adv. Neural Inf. Process. Syst, 2019, 32, 8026–8037.
48. J. Bergstra, D. Yamins and D. Cox, Proc. 30th Int. Conf. Mach. Learn., 2013, pp. 115–123.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d3cp04441f

---