Unveiling the impact of exchange-correlation functionals on the description of key electronic properties of non-fullerene acceptors in organic photovoltaics

Non-fullerene electron acceptors have emerged as promising alternatives to traditional electron-acceptors in the active layers of organic photovoltaics. This is due to their tunable energy levels, optical response in the visible light spectrum, high electron mobility, and photochemical stability. In this study, the electronic properties of two representative non-fullerene acceptors, ITIC and Y5, have been calculated within the framework of density functional theory using a range of hybrid and non-hybrid density functionals. Screened range-separated hybrid (SRSH) approaches were also tested. The results are analyzed in light of the previously reported experimental outcomes. Specifically, we have calculated the oxidation and reduction potentials, fundamental and optical gaps, the highest occupied molecular orbital and lowest unoccupied molecular orbital energies, and exciton binding energies. Additionally, we have investigated the effects of the medium dielectric constant on these properties employing a universal implicit solvent model. It was found that hybrid functionals generally perform poorly in predicting oxidation potentials, while non-hybrid functionals tend to overestimate reduction potentials. The inclusion of a large Hartree–Fock contribution to the global or long range was identified as the source of inaccuracy for many hybrid functionals in predicting both redox potentials and the fundamental and optical gaps. Corroborating with the available literature, ∼ 50% of all tested functionals predicted very small exciton binding energies, within the range of ± 0.1 eV, that become even smaller by increasing the dielectric constant of the material. Finally, the OHSE2PBE and tHCTHhyb functionals and the optimal tuning SRSH approach emerged as the best-performing methods, with good accuracy in the description of the electronic properties of interest


I. INTRODUCTION
In recent years, organic photovoltaics (OPVs) have attracted significant attention as a promising solar energy conversion device due to their low-cost, lightweight, and flexible nature. 1 In contrast to silicon-based solar cells, OPVs convert sunlight into electricity by using organic materials like polymers, small molecules, or nanocrystals as the active layer.Among such organic materials, the so-called non-fullerene electron acceptors (NFAs) have recently attracted a great deal of attention.This is mainly due to their tunable energy levels, optical absorption in the visible range, high electron mobility, and excellent photochemical stability. 2n this context, Density Functional Theory (DFT) has become a powerful tool for the computational modeling and design of novel OPV materials, providing a cost-effective and computationally efficient approach to predict electronic properties such as electron affinity (EA), ionization potential (IP), redox potentials, electronic absorption spectrum, optical gap, exciton binding energy, polarizabilities, etc. 3 However, the accuracy of DFT depends significantly on the choice of the exchange-correlation (XC) functionals.There are several such XC functionals available for electronic structure calculations, each with its own strengths and weaknesses, making the choice of the functional a nontrivial task.
In this study, we have systematically benchmarked the performance of many different XC functionals to predict some key electronic properties of two representative NFA molecules (ITIC and Y5) with an internal Acceptor-Donor-Acceptor (ADA) character (see Fig. 1).They are both non-fullerene electron acceptors that have been extensively studied for use in organic solar cells. 4ITIC stands for indacenodithieno [3,2-b]thiophene-2,8-dione, and Y5 is the commonly used label for the molecule bis(4-(dip-tolylamino)phenyl)phenylphosphine oxide.The intramolecular donor-acceptor (D-A) structures of ITIC and Y5 make them the quintessential NFAs with an internal ADA character.This kind of system can be constructed by modular synthesis, which mainly includes three elements: aromatic core (with donor character), endcapping group (with acceptor character), and alkyl side chains. 5herefore, ITIC has a multifused electron-rich core (D) flanked with two peripheral electron-deficient moieties (A), while Y5 has a core that is still electron-rich but with a DA'D character.][8][9] The properties of interest in this study are the oxidation (Φox) and reduction (Φ red ) potentials and the fundamental (ΔΦ) and optical (ΔEopt) gaps.These properties were benchmarked for ITIC and Y5 under the scope of the Density Functional Theory (DFT) and time-dependent DFT (TDDFT), using non-hybrid (pure) and hybrid generalized gradient approximation (GGA) and meta-GGA XC functionals.The acronyms GGA and meta-GGA stand for generalized gradient approximation and meta generalized gradient approximation, respectively.GGA functionals incorporate both the electronic density and its gradient in the XC term, while meta-GGA functionals also consider the kinetic energy density.Pure functionals do not include any Hartree-Fock (HF) exchange contribution, while hybrid functionals contain a certain percentage of HF exactexchange, which can be either fixed or variable based on interatomic distances, giving rise to range-separated hybrid (RSH) functionals. 1 Particularly with these functionals, the Hartree-Fock (HF) exchange rate can be finely tuned within the short and long ranges of the interelectronic density distribution.In the case of screened RSH approaches, 10,11 tuning can be employed by incorporating environmental effects into the functional itself, optimizing the accuracy of the functional in predicting the properties of the system under its dielectric environment.
To assess the accuracy of the different DFT functionals in predicting the properties of interest, we have compared the results of the DFT calculations to experimental measurements in film conditions. 12,13Additionally, we have investigated the limitations of assessing the redox potentials as simply HOMO (highest occupied molecular orbital) and LUMO (lowest unoccupied molecular orbital) energies.Furthermore, we have analyzed the impact of XC functionals on the exciton binding energies and the influence of the dielectric constant of the medium on these properties.Our results provide valuable insights into the suitability of different XC functionals in predicting the electronic properties of NFAs and could serve as a useful guide for studying and designing new OPV materials.

II. THEORETICAL MODELING
In total, we have tested 55 DFT functionals, including eight generalized gradient approximation (GGA) functionals, seven meta-generalized gradient approximation (MGGA) functionals, 26 hybrid-GGA functionals (16 short-range and 10 range-separated), and 14 hybrid-MGGA functionals (12 short-range and 2 rangeseparated).A complete list of all tested DFT functionals is shown in Table I.Additionally, we have employed a tuning protocol 11 based on the screened range-separated hybrid (SRSH) LC-wHPBE functional, as discussed in Sec.II B. The Pople basis set 6-311G(d,p) 14 was adopted for all atoms in the calculations.The geometries of ITIC and Y5 were optimized for each DFT functional in the gas phase, and harmonic vibrational analysis was performed to assure the structural stability of the optimized geometries.To save computational time, the original alkyl side chains of ITIC and Y5 were substituted with methyl groups.All calculations were performed using the Gaussian 16 program (Rev C.01). 15

A. Redox potentials and optical gap
According to the Born-Haber thermodynamic cycle approximation, [16][17][18] the Gibbs free energy variation of the Red : are given by and The different ΔG (solv) terms in Eqs. ( 3) and (4) account for the solvation free energies of the neutral (Xg), reduced (X − g ), and oxidized (X + g ) species, while the components ΔrG ox g and ΔrG red g are the free energy variations of the oxidation and reduction reactions in the gas phase, respectively, which are calculated as with the electron chemical potential, at the vacuum level, set to the zero-energy value.Moreover, the oxidation-reaction free energies in Eqs.
(3) and ( 5) are calculated in the direction of the reduction of the oxidized state, X + , to keep the commonly used convention of the redox potential calculations.The gas phase Gibbs free energy of each molecular system can be written as The first three terms on the right-hand side of Eq. ( 7) compose the internal energy of the molecule, with E, ZPE, and U 298 corresponding to the electronic total energy, zero-point energy, and thermal correction to the internal energy (including vibrational, translational, and rotational contributions), respectively.By adding the pV term to the latter ones, we obtain the molecular enthalpy.Finally, the Gibbs free energy is calculated by subtracting the term TS 298 from the enthalpy.Here, the entropic term, S 298 , also includes the vibrational, rotational, and translational contributions.
The oxidation (Φox) and reduction (Φ red ) potentials are then written in terms of the corresponding reaction free energies in solution such as where n is the number of electrons changed in the reaction and F is the Faraday constant (96.5 KJ mol −1 ).All the quantities composing the gas phase Gibbs free energy were calculated at the DFT level.The solvation free energies were obtained using the polarizable model of solvent SMD, 19 which is a "Universal Solvation Model Based on Solute Electron Density and on a Continuum Model of the Solvent Defined by the Bulk Dielectric Constant and Atomic Surface Tensions."In order to reproduce the environmental conditions of non-fullerene ITIC and Y5 thin films, the dielectric constant (ε) was obtained by the Clausius-Mossotti equation where ρ, M, N A , and α are the density of the material, molecular mass, Avogadro number, and isotropic component of the molecular polarizability, respectively. 20The ωB97XD/6-31G(d,p) theory level was chosen for this procedure due to its proven accuracy for this purpose. 21The calculated values for ITIC and Y5 resulted in dielectric constants of ε = 4.47 and ε = 4.96, respectively.These computed dielectric constants align well with the values reported in recent experimental literature for ITIC, Y6, and similar compounds, which range between 4.3 and 5.7. 21,22he fundamental gap (ΔΦ) was obtained as the difference between the oxidation and reduction potentials.The optical gap (ΔEopt) was calculated using the TDDFT. 23Here, ΔEopt is defined as the energy difference between the ground state (GS) minimum and the first Franck-Condon state (not optimized).As before, the environmental effects were accounted for via the SMD model.As defined in Eq. ( 8), the redox potentials will assume positive values with zero energy at the vacuum level.Exceptions happen when the reactionfree energies, as defined in Eqs. ( 5) and ( 6), assume positive values.In these cases, the molecule either undergoes spontaneous oxidation The Journal of Chemical Physics ARTICLE pubs.aip.org/aip/jcp or its reduction reaction is not spontaneous from the thermodynamics viewpoint.However, throughout this manuscript, we will use −ϕ i as a convention to make a consistent comparison with the calculated (through the DFT-self consistent approach) and measured HOMO and LUMO energies.These energies are usually presented as negative values.

B. Screened range-separated hybrid approach
The development of screened range-separated hybrid (SRSH) functionals has recently emerged as a promising solution to account for the impact of the electrostatic environment on the electronic properties of organic semiconductor materials. 10,11,24A general scheme for partitioning the Coulomb repulsion is the one designed for the CAM-B3LYP functional. 25It involves partitioning the Coulomb interaction using a three-parameter (α, β, ω) dependent error function for the electron interaction distance (r), where α and α + β account for the fraction of Fock exchange in the short-and long-range, respectively, and ω controls the transition from short-to long-range interactions.In the gas-phase, α + β is set to 1 to achieve the correct asymptotic potential.The parameter α is commonly set to 0.2 in order to reduce DFT selfinteraction errors, [26][27][28][29] while ω is obtained from the minimization of the following error function: where IP stands for ionization potential and EA stands for electron affinity.Such a procedure is based on the ionization-potential theorem for the neutral and anionic states of the molecular systems.When J(ω) minimization is performed in vacuum (with α = 0.2 and β = 0.8) and environment effects are subsequently incorporated through an implicit dielectric based model, polarizable continuum model (PCM), 30 this approach will be referred to as RSH-PCM, as previously proposed in the literature. 11Therefore, when J(ω) minimization is conducted in the presence of the PCM (and keeping the values α = 0.2 and β = 0.8), the method will be named RSH-PCM * .In the SRSH approach, α + β = 1/ε to maintain consistent dielectric screening with the dielectric model. 11Within this approach, the SRSH exchange-correlation functional takes the following form: 10 In the first approximation level (SRSH-PCM), 31 the relevant parameters assume the following values: ω = ω(vacuum optimization), α = 0.2 and β = 1/ε − α.In a more elaborated method (OT-SRSH-PCM), ω = ω(vacuum optimization), β is tuned based on an expression similar to that used for ω tuning, and respecting the relation α + β = 1/ε.Lastly, if new ω is obtained in vacuum using the previously optimized α, α and β are iteratively reoptimized using the latest ω (carrying out OT-SRSH-PCM) until the convergence of all parameters, the method is named OT-SRSH-PCM * (OT stands for optimal tuning).We recommend the readers to see Ref. 11 for more details about this methodological protocol.
By utilizing the tunable LC-wHPBE 32 functional combined with the basis set 6-311G(d,p) and the PCM model, we have tested all these approaches, viz., RSH-PCM, RSH-PCM * , SRSH-PCM, and OT-SRSH-PCM * .The geometry optimizations were performed at the B3LYP/6-311G(d,p) theory level, whereas the parameters of Eq. ( 9) used for the dielectric constant calculations were obtained at the ωB97X-D/6-31G(d) theory level.First, we computed the electronic energy of ionic species using vertical processes without any structural relaxation.Then, subsequently, we considered structural relaxation effects.We adopted the α, β, and ω parameters from the OT-SRSH-PCM * approach to reoptimize the geometry of neutral and ionic species and recalculate the energy levels, referred to as OT-SRSH-PCM * _R.Finally, we have incorporated the thermal corrections to compose the full Gibbs free energy calculated on the relaxed geometries to obtain the oxidation and reduction potentials.This last approach (here named OT-SRSH-PCM * _RG) is the most adequate one to make comparisons with the redox potentials measured with cyclic voltammogram experiments.

III. RESULTS AND DISCUSSION
The properties of interest in this study were the oxidation (Φox) and reduction (Φ red ) potentials and the fundamental (ΔΦ) and optical (ΔEopt) gaps of Y5 and ITIC, two representative nonfullerene acceptors.To establish a systematic comparison between theory and experiment, we selected the studies developed by Fei et al. 12 and Yuan et al. 13 as experimental references for ITIC and Y5 properties, respectively.The reported experimental values for −Φox, −Φ red , ΔΦ, and ΔEopt are −5.64,−3.92, 1.72, and 1.759 eV (705 nm) for ITIC and −5.55, −3.87, 1.68, and 1.583 eV (783 nm) for Y5, respectively.
We computed the properties of interest using various DFT functionals from different classes (hybrid and non-hybrid MGGA and GGA), as described before.The calculated values of all properties are displayed in Table I, while Figs. 2 and 3 show the differences between the calculated and experimental reference values.These differences, termed oxidation-potential (εox), reduction-potential (ε red ), fundamental (ε fund ), and optical gap (εopt) errors, are used to evaluate the accuracy of the functionals in describing the properties of interest.The smaller the difference between the calculated and experimental values, the better the XC functional's performance.First, we discuss the performance of the default functionals listed in Table I, and in Sec.III A 1, we discuss separately the performance of the tuned range separated LC-wHPBE using the protocols described in Sec.II B.
Moreover, the impact of the XC functional choice extends beyond predictions of the properties of interest.We examine (i) the relationship between redox potentials (oxidation/reduction) and HOMO/LUMO energies in light of the Koopmans' theorem; (ii) the exciton binding energy is the difference between fundamental and optical gaps; and (iii) the influence of the medium's dielectric constant on these properties, using one representative functional from each XC class.

A. Redox potentials and fundamental gap
Our analysis reveals that more than 80% of the DFT functionals that we have tested here overestimate the reduction and oxidation potentials of ITIC and Y5.In other words, the calculated values of the oxidation and reduction potentials are generally A direct correlation is observed between the performance of certain density functionals and the excessive contribution of Hartree-Fock (HF) exact-exchange in the global or long ranges, particularly in the prediction of reduction potentials.It should be clarified that the term "global" range used throughout this text means the sum of the "short" and "long" ranges.For example, several global hybrid functionals such as M06HF (100% HF), M08HF (52.2%),BHANDH (50% HF), BHANDHLYP (50% HF), M06-2X (54% HF), and M05-2X (56% HF), have a high HF contribution in the global range and are responsible for overestimated reduction potentials in this benchmarking.A similar correlation is also obtained among long-range corrected (LRC) functionals, which consider 100% HF in the long-range region, such as M11, ωB97XD, and LC-wPBE.
These functionals with high HF exchange in the global or longrange regions overestimated the reduction potential of the studied NFAs by at least 0.2 eV and up to 0.9 eV compared to experimental references.Interestingly, the inclusion of high HF exchange in the global or long-range regions does not significantly affect the oxidation potential.In fact, LRC functionals such as ωB97XD, LC-wHPBE, MN15, or global hybrid functionals such as M06-2X and M08HF can reasonably describe the oxidation potentials of ITIC and Y5 within an error margin of 0.2 eV.The results of this study indicate that a high contribution of nonlocal Hartree-Fock in a fixed amount (hybrid functionals) or in a long-range separated scheme does not benefit the redox potential calculations, especially for the reduction potentials.
By design, within the Hartree-Fock theory, the exact-exchange energy term cancels the spurious self-interaction error present in the Hartree-coulomb contribution to the electronic total energy. 33uch self-interaction errors are re-incorporated into standard functional approximations in the implementation of DFT.Therefore, the hybrid approach of combining a fraction of the HF exact-exchange energy with a conventional semi-local functional is a pathway to remedy such errors.These hybrid functionals tend to decrease unwanted electron self-repulsion and potentially alleviate the issue of over delocalized electrons often seen in LDA or GGA. 34Therefore, DFT functionals with a suitable HF contribution in the XC term can provide more accurate predictions of both oxidation and reduction potentials by incorporating both the HF contribution and electron correlation effects.
In this benchmarking, it is worth noting that the best functionals for describing reduction potential are not necessarily the best ones for describing oxidation potential.However, some functionals that give values within an error margin of 0.2 eV compared to experimental references show reasonable performance in describing both potentials.For ITIC, the functionals OHSE2PBE, M11L, and tHCTHyb are the most appropriate for describing both potentials, with associated errors ranging from −0.20 to 0.06 eV.For Y5, the functionals OHSE2PBE, tHCTH, M11L, HCTH407, tHCTHhyb, and MN12SX are the most appropriate, with errors ranging from −0.19 to 0.08 eV.
Overall, M11L, OHSE2PBE, and tHCTHhyb are the best suited for describing both the oxidation and reduction potentials of ITIC and Y5 simultaneously.The M11L is a 44-parameter local non-hybrid MGGA functional with a local character and 0% HF exact-exchange that employs dual-range local exchange to provide broad accuracy for both single-configurational (closed-shell) and multiconfigurational (open-shell) molecules and for solidstate lattice constants.The OHSE2PBE, also known as HSE03 (Heyd-Scuseria-Ernzerhof 03), is a hybrid GGA screened-exchange functional.The exchange component of the electron-electron interaction is separated into a short Perdew-Burke-Ernzerhof (PBE) and a long-range (PBE) part, and the screening parameter is fixed at a value of 0.3.Additionally, it includes 25% of HF exact-exchange in the short range.Finally, the t-HCTHhyb is a 17-parameter hybrid MGGA functional based on the t-HCTH functional, added by 15% HF exact-exchange.The t-HCTH is a non-local functional that includes the kinetic-energy density to simulate delocalized exchange effects.
The success in describing the difference between Φ red and Φox, which is known as the fundamental gap (ΔΦ), does not depend directly on the absolute values of the redox potentials.Several functionals that failed to accurately describe the redox potentials had success in describing the fundamental gaps of ITIC or Y5, as seen in Figs. 2 and 3.However, the majority of hybrid functionals tend to overestimate the experimental references of ΔΦ for both ITIC and Y5, while the majority of non-hybrid functionals tend to underestimate them.Essentially, these errors arise from hybrid functionals overestimating Φox, while non-hybrid functionals overestimate Φ red .
It is worth noting that most of the best functionals for describing ΔΦ belong to the H-GGA class, accounting for over 50% of those with errors below 0.1 eV.The H-MGGA TPSSh (ε fund = −0.07eV for ITIC and ε fund = −0.03eV for Y5), tHCTHyb (ε fund = 0.04 eV for ITIC and ε fund = 0.05 eV for Y5) and MN12SX (ε fund = 0.10 eV for ITIC and ε fund = 0.09 eV for Y5) are exceptions for describing very well the fundamental gap of both ITIC and Y5.In the top 5, the following functionals describe accurately the fundamental gap of both ITIC and Y5: OHSE1PBE (ε fund = 0.02 eV for ITIC and The Journal of Chemical Physics ARTICLE pubs.aip.org/aip/jcpε fund = 0.04 eV for Y5), HSEH1PBE (ε fund = 0.02 eV for ITIC, εfund = 0.07 eV for Y5), O3LYP (ε fund = −0.03eV for ITIC, ε fund = 0.02 eV for Y5), OHSE2PBE (ε fund = 0.04 eV for ITIC, ε fund = 0.04 eV for Y5), and tHCTHhyb (ε fund = 0.04 eV for ITIC, ε fund = 0.05 eV for Y5).Based on these results, it can be concluded that only the functional OHSE2PBE can reproduce the experimental values of Φ red , Φox, and ΔΦ for both ITIC and Y5 within an error margin of 0.1 eV.When this criterion is made a little more flexible, up to 0.2 eV, another functional stands out, the tHCTHyb.

On the performance of RSH and SRSH methods
In Fig. 4, the performance of both the RSH and SRSH methods is depicted with respect to ITIC and Y5 molecules.Notably, a significant gap renormalization is observed in the performance trends, primarily attributed to the stabilization of the electron affinity (EA) energy.When employing the vertical approximation, the functional LC-wHPBE, under its default settings (ω = 0.4, α = 0, β = 1) and using the vertical transition approach, as represented by DEF-PCM in Fig. 4, tends to overestimate the fundamental gaps by 1.29 eV for ITIC and 2.15 eV for Y5.This discrepancy arises due to the deep IP and shallow EA energy levels.It is evident that a mere optimization of the range-separation parameter (ω), as in RSH-PCM and RSH-PCM * , brings about a considerable enhancement in agreement with experimental results.However, even with these improvements, the EAs and the experimental reduction energies remain in poor agreement at these theory levels.At the OT-SRSH-PCM * theory level, the IP/EA energies of ITIC and Y5 are −5.72/−3.58eV and −5.64/−3.56eV, respectively.Therefore, besides the excellent agreement between the calculated IP and the experimental oxidation potentials, the calculated EAs are still more than 0.3 eV less negative than the experimental reduction potentials.
Then, by performing the geometry relaxation of the neutral and ionic species using the tuned LC-wHPBE functional under the OT-SRSH-PCM * approach, we obtained the energy levels represented by OT-SRSH-PCM * _R in Fig. 4. For ITIC, the geometry relaxation did not affect the IP levels, but it stabilizes the EAs by 0.11 eV.For Y5, the geometry relaxation had a smaller impact on the EA stabilization, only 0.08 eV, but changed the IP by 0.11 eV.In any case, the description of the energy levels and the fundamental gap is really benefited by the geometry relaxation.The agreement between theory and experiment at the OT-SRSH-PCM * _R level is within a negligible error margin of 0.01 eV for IP and within an error margin of 0.25 eV for EA.
It is noteworthy that the most rigorous theoretical approach for comparison with experimental measurements of redox potentials obtained through cyclic voltammetry techniques should encompass not only the electronic energy of the relaxed neutral and ionic species but also their complete Gibbs free energy.This involves consideration of zero-point energy, thermal corrections to internal energy, and contributions from molecular enthalpy and entropy, as previously delineated in Eqs.(7)  and (8).Consequently, optimal theoretical outcomes within the SRSH framework are achieved by calculating the complete Gibbs free energies associated with the oxidation and reduction reactions of the systems.The best performance is obtained at the highest level of the tuning approach, which we named OT-SRSH-PCM * _RG.
It is particularly interesting to see that these single-molecule calculation methods, further extended to include structural relaxation and thermal corrections to Gibbs free energy, are sufficient to achieve an expressive gap renormalization and excellent agreement with the experimental outcome.Such results actually place the OT-SRSH-PCM * _RG approach among the three best ones in the current study, together with the ones based on OHSE2PBE and tHCTHhyb functionals.

FIG. 4.
Energy levels of ITIC and Y5 obtained from range-separated hybrid (RSH) and screened RSH (SRSH) approaches.DEF stands for the default LC-wHPBE functional.Experimental values 12,13 of oxidation and reduction are shown for comparison.Methods' nomenclatures are defined in the methodology section.The oxidation and reduction potentials represent the energy required to remove or add an electron to a molecule, respectively.Accurate calculations of these potentials require determining the Gibbs free energy of the neutral and ionic species involved in the reactions, as it has been formulated earlier.Alternatively, the Koopmans' theorem offers a simplified approach to predict oxidation and reduction potentials by utilizing the energies of the frontier molecular orbitals obtained from the self-consistent field (SCF) convergence of the neutral molecule.According to Koopmans' theorem, the ionization potential of a molecule, which corresponds to the energy required to remove an electron from the molecule, is equal to the negative energy of the highest occupied molecular orbital (HOMO).Similarly, the reduction potential of a molecule is approximately equal to the negative energy of the lowest unoccupied molecular orbital (LUMO).It should be noted that these predictions are approximate and rely on the assumption that the electron removal or addition occurs from or to a specific orbital without any additional orbital relaxation effects.In this section, we discuss the validity of this approach in predicting the redox potentials for the studied compounds, considering all the XC functionals under investigation.

The Journal of Chemical Physics
Figure 5 shows the comparison between the calculated −Φox, −Φ red , ΔΦ and the energy of the HOMO (E HOMO ) and LUMO (E LUMO ) and their difference (E LUMO − E HOMO ) for ITIC and Y5.The orbital energies were obtained from the simple convergence of the self-consistence field (SCF) of the neutral molecules, and the redox potentials were obtained from the Born-Haber thermodynamic cycle approximation, both procedures at the same theory levels, as described before.See the supplementary material for more detailed information.We opted for suppressing the RSH and SRSH results in this discussion since both approaches are designed to obey the ionization-potential theorem, making further discussion redundant.
The purpose of this analysis is to evaluate the accuracy of approximating oxidation energy as HOMO energy and reduction energy as LUMO energy for the studied NFAs.In general, the module of the differences between Φox and E HOMO and Φ red and E LUMO for the different classes of functionals ranges from 0.02 to 0.45 eV for GGA, from 0.03 to 0.44 eV for MGGA, from 0.03 to 2.82 eV for H-GGA, and from 0.00 to 2.77 eV for H-MGGA.This indicates at first glance that the approximation of HOMO and LUMO as Φox and Φ red , respectively, is more suitable within the GGA and MGGA functionals.Among the hybrid functionals, it is notorious that those ones with a high contribution of HF in the global range [MN15 (44%), M08-HX (53.2%),M06-HF (100%), M06-2X (54%), M05-2X (56%), BMK (42%), SOGGA11-X (40.2%),BHANDH (50%), BHANDHLYP (50%)] or in the long-range (M11, ωB97X, ωB97XD, LC-wHPBE, LC-wPBE, CAM-B3LYP) are responsible for the largest disagreements within the hybrid functionals, where the redox potentials and the frontier molecular orbital energies can differ by more than 1.0 eV.If we do not consider these functionals among the hybrid functionals, the module of the differences between Φox and E HOMO and Φ red and E LUMO for the hybrid functionals ranges from 0.0 to 0.93 eV for H-GGA and from 0.0 to 0.98 eV for H-MGGA, which is still higher compared to results obtained with non-hybrid functionals.
Except for CAM-B3LYP, which uses a Coulomb-attenuating method for long range interactions, all the other long range corrected functionals tested have 100% Hartree-Fock (HF) interactions in the long range.Therefore, the inclusion of a high percentage of HF exchange-correlation in the global range or only in the long range has a considerable impact on the molecular orbital energies.In general, the HOMO energy becomes more negative than Φox, and the LUMO energy becomes less negative than Φ red , and consequently, the fundamental gap predicted as the HOMO-LUMO energy difference becomes wider than the fundamental gap predicted as the redox potential difference.For instance, for the functionals M06-HF and LC-wHPBE, the difference between such gaps is 4.75 and 4.95 eV for ITIC, respectively.
In particular, the best-performing functionals to describe the redox properties of the studied systems, tHCTHhyb and OHSE2PBE, show a reasonable agreement between the calculated redox potentials and the HOMO and LUMO energies (see arrows in Fig. 5).Specifically, using the OHSE2PBE functional, we have obtained HOMO and LUMO energies of −5.56 and −3.88 eV for Y5, and −6.30 and −3.90 eV for ITIC, respectively, which correspond to redox potentials of −5.63 and −3.91 eV for Y5, and −5.70 and −3.94 eV for ITIC.Using the tHCTHhyb functional, we have obtained HOMO and LUMO energies of −5.41 and −3.60 eV for Y5, and −6.12 and −3.63 eV for ITIC, respectively, which correspond to redox potentials of −5.41 and −3.68 eV for Y5, and −5.47 and −3.71 eV for ITIC.Therefore, the OHSE2PBE and tHCTHhyb functionals give HOMO and LUMO energies that are quite close to the oxidation and reduction potentials, respectively, for Y5.In contrast to Y5, we found that both the OHSE2PBE and tHCTHhyb functionals show a good agreement between the LUMO energy and the reduction potential for ITIC but an expressive disagreement between the HOMO energy and the oxidation potential.In fact, as shown in Fig. 5(a), all the functionals predict a HOMO energy that is more negative than the oxidation potential for ITIC, with a difference of at least 0.3 eV.This suggests that the oxidation of ITIC may be more difficult to achieve than predicted by the simple SCF energy calculations and highlights the need for caution when interpreting the redox properties of this molecule using DFT.
The ability to predict the energy levels of a system by simply calculating the frontier molecular orbital energies can be advantageous, particularly for large polymers with hundreds of atoms, where other approaches may be computationally prohibitive.However, this study highlights the importance of careful consideration when approximating redox potentials simply as the HOMO and LUMO energies, especially for functionals with high XC Hartree-Fock contributions, either in the global or only in the long ranges.Our findings suggest that such approximations may not be accurate for some functionals and that caution must be exercised when interpreting redox properties based solely on calculated HOMO and LUMO energies.

C. Optical gap
Concerning the optical gap (ΔEopt), here defined as the energy of the first excited state (S1) in the reference of the ground state (GS), Fig. 6 depicts the optical-electronic spectra of ITIC and Y5, where the most intense absorption band comes from the electronic The Journal of Chemical Physics excitation from GS to S1 state.As in both systems, this excitation is dominated by the HOMO and LUMO.In Fig. 6(b), the transition density matrix (TDM) of the GS to S1 excitation shows that this excitation presents both local centered (LC) and charge transfer (CT) contributions.One can see that a significant amount of electronic density is transferred from the NFA core [group 2 in Fig. 6(c)] to their ending groups [groups 1 and 3 in Fig. 6(c)], 45.9% in the case of ITIC and 44.7% in the case of Y5, using the OHSE2PBE method.
The local contributions to the excitation are calculated as 38.6% and 38.4% for ITIC and Y5, respectively.As with the fundamental gap, the calculated optical gaps (ΔEopt) are usually overestimated by the hybrid functionals and underestimated by the non-hybrid ones (see Figs. 2 and 3).Among the top 5, the functionals N12SX, OHSE2PBE, B3LYP, HSEH1PBE, and OHSE1PBE have the best performance in describing the ITIC optical gap.Regarding Y5, the functionals TPSSh, O3LYP, MN12L, tHC-THhyb, and M11L are at the top of the list.Describing both ITIC and Y5 within reasonable deviations from the experimental references, stand out a shorter list of functionals: tHCTHyb (εopt = −0.08 eV for ITIC and εopt = 0.06 eV for Y5), OHSE2PBE (εopt = −0.01eV for ITIC and εopt = 0.13 eV for Y5), TPSSh (εopt = −0.13eV for ITIC and εopt = 0.02 eV for Y5), O3LYP (εopt = −0.13eV for ITIC and εopt = 0.02 eV for Y5) and B3LYP (εopt = 0.01 eV for ITIC and εopt = 0.14 eV for Y5).It is important to note that 0.1 eV is equivalent to ∼50 nm in the spectral region where the low-lying absorption bands are measured (around 700-800 nm).In units of wavelength (nm), the values obtained with tHCTHyb and OHSE2PBE are red-shifted from the experimental values at 35.5 and 3.3 nm for ITIC, respectively, and blue-shifted at 28.8 and 58.9 nm for Y5, respectively.
It is, therefore, challenging to have a good description, at the same time, of both the redox potentials and the fundamental and optical gaps of such representative NFAs.In the case of redox potential assessment, we need to describe well the states, which are more ground-state related properties, while in the assessment of the optical gaps, we need to describe well the excited states.Although many functionals have been listed as good options for describing optical properties, they have almost all failed to describe the redox properties well.However, two exceptions should be noted: the HGGA OHSE2PBE and the H-MGGA tHCTHyb.Under the considered approximations adopted for the calculation of the studied properties, these functionals presented excellent performance in the description of the electronic properties of ITIC and Y5.
Furthermore, within the realm of the SRSH methods, the optical gap calculations at the OT-SRSH-PCM * _R theory level exhibit remarkable consistency.Specifically, the results demonstrate excellent alignment for ITIC (ΔEopt = 1.83 eV, in comparison to the experimental value of 1.76 eV) and reasonable concordance for Y5 (ΔEopt = 1.79 eV, compared to the experimental value of 1.58 eV).9][80] In some systems, these deformations, for instance, can have the potential to break molecular symmetry, causing a shift in the optical gap towards lower energy levels. 80In the case of Y5, addressing and refining the description of such effects through SRSH methods may be necessary to enhance the overall theoretical/experimental agreement.Furthermore, this discussion is valid in the case of any XC functional.In principle, the solvatochromic shifts from vacuum to the film dielectric constant of ITIC and Y5 are quite dependent on the XC correlation functional and the amount of HF exchange (see Tables S3 and S4).This highlights the importance of calibrating the amount of DFT and HF exchange in the short and long ranges of range-separated functionals.Exciton binding energy (E BE ) is one of the key properties that governs the physics of opto-electronic devices.The exciton binding energy is a measure of the energy required to separate the electron and hole, which have formed the exciton upon photoexcitation, into free carriers.It is commonly defined as the difference between the fundamental (ΔΦ) and optical (ΔEopt) gaps, 81,82

The Journal of Chemical Physics
Figure 7 shows the E BE calculated for ITIC and Y5 using all the tested DFT functionals.More detailed information is available in the supplementary material.Approximately 50% of all tested functionals predicted the E BE for both ITIC and Y5 within the range of ±0.1 eV, and out of this range, some functionals gave values for the E BE between −0.7 and 0.2 eV for ITIC and −0.5 and 0.3 eV for Y5.Due to the higher number of hybrid functionals available, the H-GGA and H-MGGA classes accommodate the most varied values and also the majority of the negative E BE values (more than 95%).For ITIC, five hybrid functionals are responsible for E BE < −0.2 eV: BHANDHLYP (−0.23),LC-wHPBE (−0.69 eV), ωB97X (−0.23 eV), BMK (−0.21 eV), and M06HF (−0.36 eV).In turn, for Y5, nine functionals are present with E BE lower than −0.2 eV: B1LYP (−0.43 eV), B3LYP (−0.33 eV), BHANDHLYP (−0.23 eV), CAM-B3LYP (−0.49eV), LC-wHPBE (−0.28 eV), LC-wPBE (−0.28 eV), X3LYP (−0.36 eV), ωB97XD * (−0.21 eV), ωB97X (−0.55 eV).It is noteworthy that the majority of these functionals belong to the range-separated class, with a Hartree-Fock contribution of 100% in the long range, or are global hybrid functionals with a high HF contribution in the global range.As previously discussed, these functionals result in a significant overestimation of both the fundamental and optical gaps.However, as the overestimation of the optical gap becomes more pronounced, it leads to negative values for the E BE .
Based on this benchmarking, the most successful functionals predict E BE values of 0.01 eV for both ITIC and Y5 using OHSE2PBE, 0.08 eV for ITIC and 0.09 eV for Y5 using tHCTHhyb, or −0.09 eV for ITIC and 0.0 eV for Y5 using OT-SRSH-PCM * _RG.The low E BE values predicted by these and other methods (see Fig. 7) align with recent literature. 12,13,83Considering the optical bandgap as given by the equation ΔEopt = 1240/λonset and the oxidation/reduction potentials estimated from cyclic voltammograms, the E BE of ITIC and Y5 are 0.11 and 0.30 eV, respectively, according to the experimental data provided by Fei et al. 12 and Yuan et al., 13 respectively.Moreover, if we consider λmax instead of λonset, these values decrease to −0.04 and 0.1 for ITIC and Y5, respectively.In general, λonset refers to the lowest-energy electronic transition, which often involves the promotion of an electron from the HOMO to the LUMO.This represents the minimum amount of energy required to initiate an electronic transition and initiate absorption.λmax, on the other hand, is determined by the energy of the peak with the highest intensity in absorption, which typically corresponds to an electronic transition with higher energy, "a hot electron."Furthermore, Xiaoyu Liu et al. 83 recently determined the exciton binding energy of ITIC to be 0.117 eV using temperature-dependent photoluminescence measurements.It should be noted that the exciton binding energy values reported for these molecules may vary depending on the experimental conditions and specific measurement techniques used.][84] One of the reasons for the success of NFAs could be their lower The Journal of Chemical Physics ARTICLE pubs.aip.org/aip/jcpE BE of the lower excited states, which enables efficient charge generation.The smaller the E BE , the lower the necessary driving force required for an exciton to split into the free charge carriers (electron and hole).Furthermore, free charge generation occurs via energy transfer, followed by hole transfer, which makes NFAs increasingly important for solar energy harvesting.

E. Environment effects
In this study, we compare the calculated properties of the studied compounds to the corresponding experimental data measured under film conditions.Although the real system is complex and consists of numerous NFA molecules that aggregate to form a thin film, we simplify the system by modeling it as a single NFA embedded in a continuous polarizable dielectric with a dielectric constant that is defined to reproduce the film environment conditions, as implemented in the SMD model.Consequently, in this approximation, the extension of the environmental effects is mainly determined by the dielectric constant of the medium.
The influence of the surrounding medium on the redox potentials, fundamentals, and optical gaps of ITIC and Y5 was analyzed under the different DFT functionals.We find that all the investigated properties are significantly affected by the surrounding environment.Specifically, the oxidation and reduction potentials are on average 1 eV lower and 0.7 eV higher, respectively, within the dielectric medium compared to vacuum.Furthermore, the fundamental gap is on average 1.6-1.7 eV wider, while the optical gap is 0.13 ± 0.2 eV (53 ± 23 nm) narrower in vacuum conditions.As a result, the exciton binding energy is ∼1.5-1.6 eV lower in the dielectric environment than in a vacuum.This highlights the crucial role of the dielectric constant of the medium in the generation of free charge carriers.
Increasing the dielectric constant of a material is known to weaken electron-hole interactions and decrease exciton binding energy (E BE ). 85,86The effects of varying the dielectric constant from 3.0 to 10.0 a.u. on ITIC and Y5 were investigated using four DFT functionals representing each of the non-hybrid and hybrid generalized gradient approximation (GGA) and meta-GGA (MGGA) The Journal of Chemical Physics ARTICLE pubs.aip.org/aip/jcpclasses: tHCTHyb, OHSE2PBE, M11L, and B97D3.Figure 8 shows the results obtained from this study.Assuming the lowest dielectric constant for ITIC and Y5, the E BE values are around 0.3 eV for ITIC and 0.4 eV for Y5, and it approaches zero for dielectric constants at around 5-7 a.u.for ITIC and 5-9 a.u.for Y5.These results highlight the significance of the dielectric constant of the medium in determining exciton binding energy and, ultimately, the generation of free charge carriers in organic semiconductor materials.Therefore, through molecular engineering, e.g., by incorporating polar side chains, the dielectric function could be modified to achieve spontaneous charge photogeneration and, therefore, organic devices with higher efficiency.

IV. CONCLUSIONS AND OUTLOOK
The performance of 55 XC density functionals for the calculation of key electronic properties of two representative NFAs within the DFT and TDDFT formalisms was analyzed against experimental references.This study revealed that the most accurate predictions for reduction potentials are offered by H-MGGA and H-GGA functionals, while MGGA and GGA functionals are more accurate for oxidation potentials.Among these, specific functionals stand out as the most reliable choices for predicting reduction or oxidation potentials with errors within a narrow range of 0.1 eV or even smaller.Notably, among the non-tuned functionals M11L, OHSE2PBE, tHCTH, HCTH407, and HCTH147 for ITIC and PBEPBE, M11L, tHCTH, OHSE2PBE, and HCTH147 for Y5 demonstrate superior accuracy in predicting these potentials.The study indicates a direct correlation between the performance of certain functionals and the excessive presence of HF exact-exchange in the global or long-range regions of the exchange-correlation functional.Functionals with a high HF exchange contribution in these regions tend to overestimate the reduction potential significantly.However, the oxidation potential is less affected by such contributions.Long-range corrected (LRC) functionals that incorporate HF exchange in the long-range region, as well as some global hybrid functionals, offer reasonable descriptions of oxidation potentials within a narrow error margin.The best performing density functionals were the H-GGA OHSE2PBE and the H-MGGA tHCTHyb, which were able to predict accurately both oxidation and reduction potentials and the fundamental gap of both molecules.On the other hand, in the SRSH framework, optimal results come from calculating complete Gibbs free energies for oxidation and reduction reactions, with the highest performance achieved using the OT-SRSH-PCM * _RG tuning approach.Therefore, our tests confirm the reliability of the optimal tuning screened procedure but indicate that to describe the oxidation and reduction potentials of ITIC and Y5, it is necessary to account for the molecular structure relaxation and Gibbs thermal corrections after the vertical tuning procedure.
Among the non-tuned functionals, specific ones are found to perform well in describing the optical gap for each NFA.For ITIC, the functionals tHCTHyb, OHSE2PBE, TPSSh, O3LYP, and B3LYP exhibit good accuracy, while for Y5, TPSSh, O3LYP, MN12L, tHCTHhyb, and M11L perform well.Specifically, the optical gaps obtained from OT-SRSH-PCM * _R demonstrate excellent alignment with the experiment for ITIC and reasonable concordance for Y5.It is notably challenging to simultaneously accurately describe the redox potentials, fundamental gaps, and optical gaps of NFAs.
The HGGA OHSE2PBE and H-MGGA tHCTHyb functionals and the OT-SRSH-PCM * _RG method stand out as exceptions, showing good performance for both electronic properties of ITIC and Y5.Thereafter, such methods, as well as 50% of all tested functionals, predict values of exciton binding energy lower than 0.1 eV, in good agreement with previous literature.By increasing the dielectric constant of the material, the exciton binding energy becomes even smaller, mainly because the fundamental gap is reduced.
Furthermore, approximating oxidation and reduction potentials using the energies of the HOMO and LUMO orbitals through the lens of Koopmans' theorem varies across different classes of functionals, with non-hybrid functionals demonstrating closer agreement.However, caution is advised for hybrid functionals, particularly those with high contributions of Hartree-Fock (HF) exchange-correlation, as these can lead to significant discrepancies between redox potentials and frontier molecular orbital energies.The optimal tuning screened procedure is, in this sense, an alternative once the functional is calibrated to obey the ionization-potential theorem.
Overall, our findings could aid researchers in selecting an appropriate XC functional for calculating the electronic properties of NFAs in OPVs, ultimately leading to the rational design of more efficient and cost-effective active layer materials for application in organic photovoltaic solar cells.

FIG. 1 .
FIG. 1.Chemical structures of the benchmarked systems: (a) ITIC and (b) Y5.The original alkyl side chains were substituted with methyl groups.Intramolecular acceptor (A) and donor (D) moieties in red and blue colors, respectively.

FIG. 8 .
FIG. 8. Effects of the dielectric constant on the redox potentials, optical gap, and exciton binding energy of ITIC (a) and (b) and Y5 (c) and (d).

TABLE I .
The oxidation (Φox) and reduction (Φ red ) potentials and the fundamental (ΔΦ) and optical (E opt ) gaps of ITIC and Y5 calculated with different DFT functionals.Energy values in units of eV.Experimental values presented at the bottom.
FIG. 2. The differences between the calculated and experimental Ref. 12 values of the oxidation and reduction potentials and the fundamental and optical gaps of ITIC.The shaded region demarcates the range between −0.1 and 0.1 eV.FIG. 3. The differences between the calculated and experimental Ref. 13 values of the oxidation and reduction potentials and the fundamental and optical gaps of Y5.The shaded region demarcates the range between −0.1 and 0.1 eV.