Glass transition temperature prediction of disordered molecular solids

Amorphous molecular solids, characterized by disordered molecular packing, are vital in various electronic applications due to their flexibility, processability, and tunability. A key property influencing the stability and lifetime of these materials in devices like organic light-emitting diodes (OLEDs) is the glass transition temperature (*T_g*). Identifying high-*T_g* compounds is crucial. While previous attempts at *T_g* prediction have employed statistical approaches like quantitative structure-property relationships (QSPR) and machine learning, as well as direct extraction from molecular dynamics (MD) simulations, a reliable and systematic approach for a diverse dataset of conjugated small molecules has been lacking. Challenges include the computational cost of forcefield parameterization and MD simulations, and the difficulty in creating a protocol adaptable to molecules with varied building blocks. The authors' previous work predicted *T_g* values significantly deviating from experimental results, indicating a need to adjust non-bonded forcefield parameters. This paper introduces a novel methodology addressing these challenges.

The paper reviews existing methods for predicting glass transition temperature (*T_g*), highlighting the limitations of previous approaches. These include QSPR and machine learning techniques, as well as methods relying on direct extraction from molecular dynamics (MD) simulations. The authors point out the computational cost and human effort involved in existing force field parameterization protocols, particularly the difficulty of adapting them to diverse molecular structures. Previous work by the authors using a similar protocol, but with force fields from the OPLS database, resulted in inaccurate *T_g* predictions, motivating the current study's focus on improved force field parameterization using non-bonded terms.

The proposed methodology consists of two parts: (1) a fitting protocol for density-temperature (*ρ*-*T*) plots and (2) non-bonded forcefield parameterization using atoms-in-molecule electron density partitioning (DDEC6). The conventional bilinear fit for *T_g* determination from *ρ*-*T* plots suffers from ambiguity in selecting fitting ranges, introducing human bias. To mitigate this, the authors introduce an *R²*-*T* plot to identify optimal fitting ranges, which provide unambiguous ranges for linear regressions. The non-bonded parameters (atomic partial charges and Lennard-Jones parameters) were derived using the DDEC6 method. The electron density was obtained using Gaussian 16 at the wB97X-D3/6-311G(d,p) level. DDEC6 computations were performed using Chargemol. The molecules were partitioned into rigid fragments, and dihedral potentials were parameterized using constrained optimization scanning performed at the wB97X-D3/def2-TZVP level using ORCA 4.2.1. Classical MD simulations were conducted using GROMACS, employing PME for long-range electrostatic interactions and a cutoff for non-bonded interactions. The system was initially generated using Packmol, heated, equilibrated, and then cooled at a rate of 100 K ns^-1. The *T_g* was then determined using the developed fitting protocol. The effect of various parameters (cooling rate, system size, etc.) was also investigated.

The study demonstrated that the proposed fitting procedure significantly reduces the human analysis variation associated with conventional bilinear fitting, decreasing the variation in *T_g* prediction from around 90 °C (for BCP) to a much smaller range. Using the DDEC6 force field parameters, the method achieved a mean absolute error (MAE) of 20.5 °C for 24 compounds (excluding two outliers), representing a substantial improvement compared to the 59.1 °C MAE obtained using the OPLS-MK forcefield. The *R²* value also increased from 0.61 to 0.87. Comparison with the protocol by Patrone et al., which also aims to reduce human bias, showed that the proposed method outperformed it in terms of both MAE (20.5 °C vs 64.7 °C) and *R²* (0.87 vs 0.75). The computational cost analysis revealed that MD simulations (compressing, heating, and cooling) consumed around 85% of the total computational resources. The DDEC6 force-field parameters were made publicly available.

The improved accuracy of *T_g* prediction using the proposed methodology highlights the importance of accurate non-bonded forcefield parameters. The significant reduction in MAE from 59.1 °C to 20.5 °C demonstrates that the combination of the improved fitting protocol and DDEC6 forcefield parameterization leads to a reliable and efficient approach for *T_g* prediction. The remaining discrepancy for the two outlier compounds suggests that improved force field models, potentially incorporating non-local correlations in dihedral potentials, are needed for molecules with specific structural features like the TPA moiety commonly found in organic semiconductors. The readily available parameters can facilitate multiscale simulations of complex systems in OLEDs and in silico deposition protocols.

This study presents a reliable and automated computational methodology for *T_g* prediction in disordered molecular solids. The combination of the improved fitting protocol and DDEC6 forcefield parameterization significantly enhances accuracy and reduces human bias. While further forcefield development is needed to address remaining discrepancies, this approach is a valuable tool for prescreening organic semiconductor materials for OLED applications. The open-access availability of the generated force field parameters will further contribute to the advancement of the field.

The study identified two outliers (MTDATA and 2-TNATA), starburst TPA-based dendrimers, where the predicted *T_g* values differed significantly from experimental data. This indicates limitations in the current forcefield's ability to accurately capture the behavior of molecules with this specific structural motif. The high computational cost of the MD simulations, particularly the heating and cooling steps, also presents a challenge for large-scale screening applications. Future work should focus on improving the efficiency of the MD simulations and developing forcefield parameters that can capture the complexities of molecules such as TPA-based dendrimers.