Singlet-Fission-Induced Enhancement of Third-Order Nonlinear Optical Properties of Pentacene Dimers

Using quantum chemical calculations and exciton dynamics simulations, we investigate the static second hyperpolarizability γ [the third-order nonlinear optical (NLO) property at the molecular scale] of slip-stacked pentacene dimer models in the correlated-triplet-pair [1(TT)] state created from the singlet excited state in the singlet fission (SF) process. It is found that the SF induces significant (∼20 times at maximum) enhancement of γ/monomer in the 1(TT) state as compared to that in the singlet ground state. The origin of the remarkable enhancement of γ/monomer is revealed by analyzing the γ density distribution and the intermolecular orbital interaction. Furthermore, we clarify molecular packings suitable for highly efficient SF and largely enhanced γ values of a novel class of SF-induced NLO systems, which have promising potential to surpass the conventional NLO systems.


INTRODUCTION
Third-order nonlinear optical (NLO) phenomena have attracted a great deal of attention in a variety of fields including physics, chemistry, biology, and materials science due to their wide applications including extremely large data storage, 1,2 ultrafast optical switching, 3 ultrahigh sensitivity spectroscopy, 4 and nanofabrication. 5,6 Because the materials with both the large third-order NLO property and ultrafast response time are indispensable for realizing such applications, many experimental and theoretical studies have been conducted to explore the mechanism of NLO properties in molecular systems as well as to construct design guidelines for highly active NLO materials. So far, several structure−NLO property relationships and design guidelines have been revealed; for example, introduction of donor/acceptor groups, 4,7 extension of the π-conjugation length, 8 and tuning the charge states 9 have been applied to the closed-shell molecular systems to enhance the second hyperpolarizability γ, which is the third-order NLO property at the molecular scale. In this regard, in recent years, open-shell singlet systems have been explored from the viewpoint of highly active NLO substances since the theoretical prediction of significant enhancement of γ for open-shell singlet molecules with intermediate diradical character y by Nakano et al. 10−12 The diradical character is a well-defined chemical index for the electronic state of open-shell singlet systems and takes a value between 0 (closed-shell) and 1 (pure open-shell). 13,14 The theoretical prediction of the y−γ correlation and the novel design principle for highly efficient open-shell NLO systems have been substantiated by accurate ab initio molecular orbital and density functional theory (DFT) calculations on model and realistic molecular systems 12,15−20 and experimental measurements including two-photon absorption (TPA) 21,22 and third-harmonic generation 23,24 for several synthesized open-shell singlet molecules.
Furthermore, Nakano and co-workers have clarified spinmultiplicity dependences of γ in molecules with different diradical characters in the singlet state, that is, the relationship between spin-multiplicity, diradical character, and γ: the closed-shell or open-shell molecules with weak diradical character exhibit the enhancement of γ with the increasing spin-multiplicity, whereas the molecules with intermediate and large diradical characters exhibit the reduction of γ with the increasing spin-multiplicity. 25 This finding suggests the possibility of another novel class of open-shell NLO materials using high spin states, which are far superior in magnitude to the conventional closed-shell NLO systems and even openshell singlet NLO systems with intermediate diradical character.
However, few realistic NLO systems utilizing high spin states have been proposed because it is known to be difficult to convert a singlet into a higher spin state due to small spin− orbit coupling in organic molecules. To overcome this problem, we focus on the singlet fission (SF) phenomenon, which is an ultrafast process creating two triplet excitons from one singlet exciton 26−31 through a spin-allowed process to a correlated-triplet-pair (singlet overall) state. The SF is now studied intensively due to its possibility to improve the photoelectric conversion efficiency in organic solar cells. 32 As the feasibility conditions of SF systems, Michl and co-workers have proposed the energy level matching conditions at the single molecular level: (i) 2E(T 1 ) ∼ E(S 1 ) or 2E(T 1 ) < E(S 1 ) and (ii) 2E(T 1 ) < E(T 2 ), where S 1 , T 1 , and T 2 indicate the lowest singlet and triplet excited states and the second lowest triplet excited state, respectively. 28,29 On the basis of these conditions, Minami and Nakano have constructed a novel design guideline for SF molecules using diradical (y 0 ) and tetraradical (y 1 ) characters: the molecules with small diradical (∼0.1 < y 0 < ∼0.5) and much smaller tetraradical character (y 1 /y 0 < ∼0.2) tend to satisfy these energy level matching conditions. 33 Thus, SF is predicted to be a feasible way to efficiently and immediately generate the triplet state as the high spin state enhancing the NLO property. Moreover, from the viewpoint of "weak diradical character", we expect that SF molecules are good candidates for a novel class of open-shell NLO systems, referred to as SF-induced-NLO systems, by combining the SF and high-spin NLO guidelines mentioned above. Such molecules have the possibility to exhibit NLO properties that surpass those of open-shell singlet molecules with intermediate diradical character. From this hypothesis, it is predicted that the SF systems have the potential for the ultrafast NLO switching systems using the drastic change in the NLO property through SF. Indeed, larger NLO responses have been reported recently in the SF process of pentacene-derivative crystals, though the NLO enhancement mechanism remains veiled. 34,35 In this study, therefore, we investigate the effects of SF on the third-order NLO properties (γ) of slip-stacked pentacene dimer models with various intermonomer configurations to clarify the guiding principle for controlling the γ in the correlated-triplet-pair [ 1 (TT)] state created by SF. The obtained results will contribute to pioneering a novel class of highly efficient NLO materials, that is, SF-induced NLO materials.

METHODOLOGY
2.1. Pentacene Dimer Models. The geometry of pentacene monomers was optimized by the spin-flip timedependent density functional theory (SF-TD-DFT) method with the Tamm−Dancoff and noncolinear approximations using the 6-311G* basis set and BHHLYP xc-functional. 36−38 We construct slip-stacked dimer models using the optimized monomer structure, where z and x axes indicate the longitudinal and lateral directions, respectively, and the y axis represents the stacking direction (the intermonomer distance is fixed to be 3.50 Å) ( Figure 1). The center of the red-colored molecule is located at the origin of the coordinate axis.

2.2.
Representation of the Electronic Structure of the 1 (TT) State Using the Spin-Unrestricted DFT Calculation. In this study, the broken-symmetry (spin-unrestricted) DFT [BS-DFT (UDFT)] calculations are employed to describe the electronic structure of the 1 (TT) state. We first calculate the triplet states for constituent pentacene monomers using the UDFT method and prepare the initial guess for the calculation of pentacene dimers, which is constructed from up and down spin wave functions localized at each monomer. Next, the electronic structure of the 1 (TT) state is obtained by the UDFT calculation with the prepared initial guess. Note here that the 1 (TT) state is expressed by the multiconfigurational determinants, which implies that the spin-restricted TDDFT method cannot describe the electronic state correctly. On the other hand, the spin-unrestricted single determinant wave function such as UHF and UDFT is known to have the ability to expand in the form of a limited configuration interaction, 39,40 though including spin contamination parts.
Indeed, the open-shell singlet ground states for various polycyclic aromatic hydrocarbons are found to be well described by UDFT methods, 15,41−43 and for their NLO properties, we have already justified the UDFT results, e.g., LC-UBLYP results, for the hyperpolarizabilities by comparing with strongly correlated results such as UCCSD(T) results. 18 In contrast, since the reliability of the UDFT 1 (TT) solutions for energies and NLO properties has not been elucidated well, we compare the total energies, excitation energies, and γ values calculated by the UDFT and the restricted active space with double spin-flip (RAS-2SF) methods, the latter of which more correctly describes the 1 (TT) state. 44−48 As a result, the variations in the γ values and the excitation energies for a geometric change at the UDFT level of approximation are found to well reproduce those at the RAS-2SF level of approximation (see Figure S4−S8 in the Supporting Information). Therefore, we discuss those quantities in the 1 (TT) state at the UDFT level of approximation in this study.
2.3. Calculation and Analysis Methods for Second Hyperpolarizability. We employ the finite-field (FF) approach 49 to evaluate the longitudinal component of the static γ (γ zzzz ). The γ zzzz value is calculated using the fourthorder numerical differentiation formula of energy E with respect to the longitudinal electric field F z where E(F z ) indicates the total energy in the presence of F z . Also, we employ the LC-UBLYP(μ = 0.33 bohr −1 )/6-311+G* method 50 for the FF approach. In the present calculations, we choose F z = 0.0005 au, which is found to achieve a numerical accuracy within an error of ∼1% on the static γ value (see Section 5 in the Supporting Information). To clarify the spatial contribution of electrons to γ zzzz , we perform γ density analysis. 51,52 The longitudinal γ density ρ zzz (3) (r) is defined as the third-order derivative of the electron density ρ(r) with respect to F z

ACS Omega
where r z represents the longitudinal component z of the electron coordinate. Also, to clarify the origin of the enhancement of γ based on the summation-over-state (SOS) expression of γ, we examine the excitation energies (E ex ) and the longitudinal transition moments (μ z ), which contribute to the γ values, using the TD-LC-UBLYP (μ = 0.33 bohr −1 )/6-311+G* method. The SOS expression of γ for symmetric systems is given by where E 0n denotes the excitation energy to the nth excited state from the initial state 0 and μ nm is the transition moment along the longitudinal (z) axis between the nth excited state and the mth excited state. 53

SF Dynamics Simulation with the Quantum Master Equation Approach
. To investigate the time evolution of 1 (TT) population [P( 1 (TT))] in each dimer configuration, we conduct SF exciton dynamics simulations with the second-order time-convolutionless (TCL) quantum master equation (QME) approach. 54−56 We consider the total Hamiltonian of the pentacene dimer models given by H = H ex + H ph + H ex-ph , where H ex is the exciton Hamiltonian; H ph and H ex-ph represent the phonon (vibration) and exciton−phonon (vibronic) coupling Hamiltonian, respectively (see Section 1 in the Supporting Information). We assume only the Holstein coupling, which is predicted to provide significant effects on the SF dynamics. 57−59 The second-order TCL QME is expressed by 43 represents the adiabatic exciton state (the eigenstate of H ex with the eigenenergy of E α ) given by a linear combination of diabatic exciton state {|m⟩} = {|S 1 S 0 ⟩, |S 0 S 1 ⟩, | CA⟩, |AC⟩, |TT⟩}. The second term on the right-hand side of eq 4 indicates the relaxation of exciton density matrix ρ ex (t), which is characterized by the relaxation rate γ m (ω,t) expressed by (under Markov approximation) 45 where J m (ν) is the spectral density of the Holstein vibrational mode of the mth diabatic exciton state and n(ν, T) indicates the Bose−Einstein distribution at the temperature T. J m (ν) is approximated by an Ohmic function with a Lorentz−Drude cutoff 56 where λ m is the reorganization energy and Ω m is the cutoff frequency in the mth diabatic exciton state. In this study, the diagonal and off-diagonal elements of the exciton Hamiltonian H ex are calculated at the LC-RBLYP/6-311G* level of approximation using the range-separating parameter μ (0.20 bohr −1 ) determined by the IP tuning method. 60,61 The excitation energies of the Frenkel exciton (FE) state and the 1 (TT) state are approximated using those at the monomer level: E(FE) ∼ E(S 1 ) = 2260 meV and E( 1 (TT)) ∼ 2E(T 1 ) = 1720 meV. The charge-transfer (CT) state excitation energy E(CT) and the electronic coupling are estimated at each dimer configuration (see Section 1 in the Supporting Information). The vibronic coupling parameters are set to (λ m , Ω m ) = (λ, Ω) = (50, 180) meV, which are known to be the typical carbon−carbon stretching mode for acene molecules. 58 The time evolution of P( 1 (TT)) is obtained by applying the six-order Runge−Kutta method to solving the QME (eq 5) with the initial state S 1 S 0 at T = 300 K (where S 0 is the singlet ground state), and the P( 1 (TT)) at 1 ps is regarded as the 1 (TT) yield in this study.
2.5. Evaluation of Energy Differences between Excited States in the Singlet Fission Process. For evaluation of the lifetime of the 1 (TT) state, we examine the energy difference between adiabatic FE and 1 (TT) states (ΔE(FE − 1 (TT)) = E(FE) − E( 1 (TT))) and that between 1 (TT) and 5 (TT) states (ΔE(Q − S) = E( 5 (TT)) − E( 1 (TT))), both of which are predicted to be related to the lifetime of 1 (TT). ΔE(FE − 1 (TT)) is estimated from diagonalization of H ex , and ΔE(Q − S) is calculated by the CAM-UB3LYP/6-311G* method, 62 which is found to well reproduce the tendency of ΔE(Q − S) calculated using the restricted active space with double spin-flip (RAS-2SF) method (see Figure S6 in the Supporting Information). 44 Moreover, using the CAM-UB3LYP/6-311G* method, we calculate the multiple diradical character y i , 52,63 which is defined as the occupation number of the lowest unoccupied natural orbital (LUNO) + i (n LUNO+i ), to evaluate the triplet− triplet interaction in the 1 (TT) state. We employ y 0 and y 1 to characterize the tetraradical 1 (TT) state, together with the average diradical character y av , which is defined by the arithmetic average of y i . 64 All of the calculations, except for the RAS-SF calculations, were performed by the Gaussian 09 program package. 65 The RAS-SF calculations were performed by the Q-Chem 4 program package. 66

RESULTS AND DISCUSSION
3.1. Spin State Dependence of γ of Pentacene Monomers. We discuss the spin-multiplicity dependence of the calculated γ values of the monomer. It is found that the γ value increases from γ(S 0 ) = 1.39 × 10 5 au to γ(T 1 ) = 5.71 × 10 5 au by changing from the singlet to the triplet state. We compare the excitation having the largest |μ z | value from the singlet ground state (S 0 → S 5 excitation) with that from the first triplet state (T 1 → T 4 excitation) (see Tables S1 and S2 in the Supporting Information). The enhancement of γ in the T 1 state is predicted to be caused by a significant reduction in E ex with keeping a large amplitude of |μ z | for the T 1 → T 4 excitation. This tendency is found to accord with the highspin open-shell NLO design guideline. 25 3.2. Molecular Packing Dependences of γ at the 1 (TT) State, 1 (TT) Yield, and SF-Induced NLO Efficiency. We examine γ at the 1 (TT) state and P( 1 (TT)) in each dimer configuration. Figure 2a,b shows the molecular configuration
3.3. Analysis of γ at the 1 (TT) State Using γ Density and Excitation Properties. To clarify the spatial contribution of electrons to the enhancement of γ, we investigate the γ density distributions for the dimer models A, B, and C ( Figure  3). It is found that the γ density distributions for A and C exhibit well-separated positive and negative densities with significant amplitudes between the monomers, which largely enhance the γ amplitudes due to the third-order field-induced intermolecular charge transfer (CT), whereas that for B does not exhibit such features, resulting in the smaller γ amplitude.
We further examine the excitation energy E ex and the longitudinal transition moment μ z of 1 (TT 1 ) → 1 (TT n ) excitations for A, B, and C to clarify the origin of the enhancement of γ. Here, 1 (TT n ) represents the nth 1 (TT) state, where the lowest 1 (TT 1 ) is also referred to as 1 (TT) as mentioned above (see Tables S4−S6 in the Supporting  Information). In general, it is difficult to obtain quantitative γ values by the SOS models in the realistic systems because it requires precise excitation energies and transition moments over a large number of excited states. On the other hand, a small number of low-lying excited states with moderate transition moments from the 1 (TT) state are expected at least to primarily contribute to the increase in γ based on the SOS formula eq 4. First, we discuss the characteristics of 1 (TT 8 ) (for C) and 1 (TT 9 ) (for A and B) states with the first largest transition moment amplitudes from 1 (TT). The primary excited configurations in these excitations are found to be the same as those in the T 1 → T 4 excitation in the pentacene monomer (see Sections 3, 5, and 6 in the Supporting Information). The highest occupied molecular orbital (HOMO) −2 → lowest unoccupied molecular orbital (LUMO) and HOMO → LUMO+2 excitations in the 1 (TT) state are found to correspond to the HOMO−1 → HOMO and LUMO → LUMO+1 excitations in the T 1 state, respectively (see, for example, Figure S10 for A in the Supporting Information). As shown in Table 1, these excited states have similar E ex and |μ z | values at A, B, and C molecular configurations. This result indicates that 1 (TT 8 ) and 1 (TT 9 ) states give similar contributions to the γ values for A, B, and C.
Next, we discuss 1 (TT 1 ) → 1 (TT 3 ) (for A and B) and 1 (TT 1 ) → 1 (TT 4 ) (for C) excitations with the second largest |μ z | values. Although these excited states are found to be dominated by the HOMO → LUMO transition, the |μ z | values

ACS Omega
Article are shown to be significantly different from each other ( Table  2). To clarify the origin of this difference, we investigate the transition density (ρ hl ) distribution between the α-HOMO and α-LUMO in the 1 (TT) state ( Figure 4). As seen from the ρ hl distribution, the HOMO → LUMO transitions for A and C correspond to the intermolecular CT with large transition density amplitudes, which lead to the large |μ z | for A and C. Moreover, this intermolecular CT transition is predicted to be the excitation from the 1 (TT) state to the CT state because the primary distributions of α-HOMO and α-LUMO have the distribution of the LUMO and HOMO in the S 0 state on the upper and lower monomers, respectively (Figures S10 and S12). In contrast, the HOMO → LUMO transition for B is found to have a small amplitude of |μ z | because the transition density is relatively localized to a small spatial region with small amplitudes, and the relative phase difference of the ρ hl in the monomers causes some cancellations in the contributions to μ z . As seen from the summation-over-state expression of γ in the perturbation theory, these larger amplitudes of the intermolecular transition moment for A and C than for B are predicted to contribute to more enhancement of γ for A and C than for B. 53 Therefore, the relative γ values for the dimers on the (z, x) plane are predicted to be determined by not the excitation with the first largest |μ z | value but that with the second largest |μ z | value, which corresponds to the HOMO → LUMO transition with a smaller orbital energy gap.
Here, we consider the electronic couplings V hl and V lh , which are defined as V ij = ⟨i|F|j⟩(F: the Fock operator) and are included in the off-diagonal elements in H ex , because the degree of the intermolecular HOMO−LUMO overlap is related to the relative values of V hl and V lh . V lh indicates the electronic coupling between the LUMO of the upper monomer (l A ) and the HOMO of the lower monomer (h B ), that is, between CT and 1 (TT) states in the SF process. Therefore, as seen in Figure S1b,c in the Supporting Information, V hl and V lh strongly depend on the molecular configurations because the electronic coupling is originated from the orbital overlap at a pair of neighboring chromophores. Since V lh is approximately regarded as the transfer integral between the l A and h B for the dimer models, we predict that if the |V lh | is large, the distributions of the HOMO and LUMO in 1 (TT) state extend over the neighboring monomer by including both components l A and h B . It is found that |V lh | exhibits large values for A and C (Table 2) because of the large in-phase or out-of-phase overlap between the l A and h B . Therefore, it is found that the α-HOMO and α-LUMO, which have primary distributions on the upper (l A ) and lower (h B ) monomers, respectively, also have slight distributions on the lower (h B ) and upper (l A ) monomers, respectively, for A and C (see Figure 5), resulting in the large positive and negative ρ hl densities well-separated on the upper (lower) and lower (upper) monomers for A and C (see Figure 4). In contrast, |V lh | shows a smaller value for B since the contributions of the in-phase or out-of-phase overlap between the l A and h B are found to be almost canceled with each other. This small |V lh | value is predicted to cause the slightly extended α-HOMO and α-LUMO distributions on the lower and upper monomers, respectively, in the 1 (TT) state for B ( Figure 5), resulting in the ρ hl limited to a small spatial region (Figure 4). For β-HOMO and β-LUMO, the same discussion can be made by replacing the upper and lower monomers.
From these results, we can see that the molecular configuration dependence of |V lh | values correlates with that of the ρ hl distribution patterns. In brief, we found the relationship between ρ hl and |V lh |, that is, the ρ hl distribution exhibits the character of the intermolecular CT transition from 1 (TT) to CT states on the dimer configurations with the large amplitude of the electronic coupling V lh between CT and 1 (TT) states in the SF exciton Hamiltonian. Now, we predict that the large |V lh | region well corresponds to the large γ region since the ρ hl distribution with significant amplitudes over the dimer and well-separated positive and negative distributions between the monomers cause the enhancement of the HOMO → LUMO transition moment amplitude. To confirm this prediction, we compare the two-dimensional maps [(z, x)− dependences] of γ and |V lh | (Figures 2a and 6). It is found that the molecular configurations with the relatively large |V hl | values tend to exhibit the large γ values, though the peaks of the V lh map do not necessarily coincide with those of the γ map. Consequently, the magnitudes of the SF electronic

ACS Omega
Article coupling V hl and V lh between CT and 1 (TT) states are found to be an effective index determining the relative amplitude of γ in the 1 (TT) state. Indeed, we evaluate the γ in the 1 (TT) state for two realistic SF dimer models (1 and 2 in Figure S15a in the Supporting Information) in a pentacene crystal. The dimer model 2 having the large displacement along the longitudinal (z) axis exhibits a larger γ value (9.33 × 10 5 au) than 1 (5.22 × 10 5 au) due to the third-order intermolecular CT (see Figure S16a in the Supporting Information). Namely, such monomer configuration control for enhancing γ is found to work well not only in the artificial slip-stacked dimer model but also in the realistic pentacene crystal with a herringbone packing. Moreover, as a possible approach to achieving the molecular configuration with large γ value and 1 (TT) yield, we could consider a chemical modification such as an introduction of appropriate substituent groups to zigzag edges of pentacene. 6,13-Bis(triisopropylsilylethynyl) pentacene (TIPS-pentacene) is known to undergo the SF and slip-stacked dimer configuration with large displacement along the longitudinal axis (see Figure  S15b in the Supporting Information). The γ value in the 1 (TT) state of the TIPS-pentacene dimer model is found to be larger (γ = 25.3 × 10 5 au) than those of the herringbone pentacene dimer models (e.g., γ = 9.33 × 10 5 au for the dimer model 2). Therefore, among realistic pentacene frameworks, we expect that the TIPS-pentacene is a good candidate for highly efficient SF-induced NLO materials. We discuss the aggregation effect in molecular aggregates on γ, though it is not the main purpose of this study. The important factor in the enhancement of γ in the 1 (TT) state is whether the HOMO and LUMO are delocalized over the neighboring monomer. As we mentioned above, the large electronic coupling leads to the delocalization of the HOMO and LUMO in the 1 (TT) state. In molecular aggregates with large electronic couplings between monomers, the HOMO and LUMO distributions are predicted to be more delocalized than in the dimer models. Since the open-shell singlet multiradicaloids with intermediate open-shell characters are predicted to significantly enhance the γ values along the intermolecular interaction pathway, 19,64,67 the γ values are expected to increase more in larger-size molecular aggregates, though the size and structural dependences of such effects should be clarified in detail.
3.4. Energy Differences between Frenkel Exciton and 1 (TT) States and between 1 (TT) and 5 (TT) States. Furthermore, from the viewpoint of realizing the control of the NLO property in the 1 (TT) state created by SF, we investigate the lifetime of the 1 (TT) state. If the lifetime is short, it is difficult to utilize the remarkably large γ in the 1 (TT) state since the 1 (TT) state is predicted to rapidly dissociate into two isolated triplet (2 × T 1 ) states via the 5 (TT) state or to recreate the FE state by triplet−triplet annihilation (TTA) (Figure 7a). Therefore, large ΔE(FE − 1 (TT)) (= E(FE) − E( 1 (TT))) and ΔE(Q − S) (= E( 5 (TT)) − E( 1 (TT))) values are expected to be one of the conditions for utilizing the 1 (TT) state. We calculate ΔE(FE − 1 (TT)) and ΔE(Q − S) at the molecular configurations with fixed z = 7.0 Å as a function of displacement x, where a γ peak is included. Figure 7b,c shows the variations in ΔE(FE − 1 (TT)) and ΔE(Q − S) as a function of x. It is found that the dissociation process more easily proceeds than the TTA process because ΔE(Q − S) is shown to be much smaller (more than 1 order) than ΔE(FE − 1 (TT)). Interestingly, ΔE(Q − S) shows the peak for C, where the γ is significantly enhanced, that is, the 1 (TT) state for C has a longer lifetime than those in the other intermonomer configurations. We examine the reason why ΔE(Q − S) has the peak for C using the average diradical character y av (= (y 0 + y 1 )/2) for 1 (TT). The stabilization of the 1 (TT) state by the intermolecular bonding interaction between the two triplet excitons is predicted to be described by the y av value since the large (small) y av value indicates the weak (strong) covalent-like intermolecular interaction. Figure 7c shows the molecular packing dependences of ΔE(Q − S) and y av . It is found that the increase and decrease in ΔE(Q − S) correspond to the decrease and increase in y av value, respectively.
We examine the x-dependence of y av using natural orbitals (NOs) corresponding to y 0 and y 1 for C and (z, x) = (7.0, 0.0) Å (D). For D, the highest occupied natural orbital (HONO) and LUNO in the 1 (TT) state are found to be made by the out-of-phase and in-phase mixings of the HOMO in the monomer S 0 state (corresponding to singly occupied natural orbital (SONO1) in the T 1 state), respectively (Figure 8a,b). Similarly, the HONO−1 and LUNO+1 are found to be made by the out-of-phase and in-phase mixings of the LUMO in the monomer S 0 state (corresponding to SONO2 in the T 1 state), respectively. This implies that there exist only the interactions between the HOMO−HOMO and LUMO−LUMO in D. On the other hand, for C, the HONO and LUNO in the 1 (TT) state exhibit primary distributions at the upper (lower) zigzag edge in the upper (lower) monomer, whereas the HONO−1 and LUNO+1 do primary distributions at the lower (upper) zigzag edge in the upper (lower) monomer (Figure 8d). Here, we consider the linear combinations of the SONO1 and

ACS Omega
Article SONO2 in the T 1 state, which are defined as the SONO1 − SONO2 and SONO1 + SONO2 (Figure 8c). In Figure 8c,d, it is found that the HONO−1 − LUNO+1 is described by the linear combinations of the SONO1 and SONO2 in the T 1 state, that is, the HOMO and LUMO in the S 0 state. This result indicates that there exist not only the HOMO−HOMO and LUMO−LUMO intermolecular interactions but also the HOMO−LUMO intermolecular interaction in C. Thus, it is found that the intermolecular interactions for C are stronger than those for D because the HOMO−LUMO intermolecular interaction does not exist in D.
The strong intermolecular interaction in C is predicted to cause the stabilization of the 1 (TT) state, resulting in the smaller y av and the larger ΔE(Q − S) values than D. For larger x values (>1.5 Å) in Figure 7c, the y av and ΔE(Q − S) exhibit increase and decrease behaviors, respectively, because of the decreasing intermolecular structural overlap. Judging from the fact that the HOMO−LUMO intermolecular interaction relates to the electronic coupling V hl and V lh , we can expect that the intermonomer configurations with large |V hl | values tend to have large ΔE(Q − S) values. From these results, the relationship among the SF efficiency, third-order NLO property, and the lifetime of the 1 (TT) state is found to exhibit the positive correlation in the molecular packing dependence, that is, the molecular packing suitable for efficient SF also tends to exhibit a large enhancement of the third-order NLO property in the 1 (TT) state with a longer lifetime. Note here that the efficient SF implies a rapid creation of 1 (TT) in this study, not implying the consecutive rapid dissociation into two isolated triplet excitons, though the SF process usually includes the dissociation process of 1 (TT) to T + T. 28,29 Recently, the separation and the lifetime of the 1 (TT) state have been experimentally measured for several SF systems. For example, it has been reported that the separation of the 1 (TT) state in rubrene single crystals occurs on the picosecond time scale. 68 Moreover, it has been found that quinoidal bithiophene (QOT2) has the long-lived 1 (TT) state, where the lifetime of the 1 (TT) state and ΔE(Q − S) are microsecond time scale and 230 meV, respectively. 69,70 From the viewpoint of using the 1 (TT) state, we expect that SF systems having the long-lived 1 (TT) state such as QOT2 can be good candidates for SF-induced NLO systems.

CONCLUSIONS
In this study, we have explored the effective molecular packing in the pentacene slip-stacked dimer models suitable for high SF efficiency [rapid creation of 1 (TT)] and the large third-order NLO property in the 1 (TT) state. For the intermonomer packing with a large |V hl | value, which relates to the 1 (TT) creation rate, a remarkable increase in the γ is predicted, the feature of which is found to be caused by the HOMO−LUMO transition in the 1 (TT) state with a small excitation energy and a large intermonomer transition moment amplitude. The γ value is found to be maximized in (z, x) = (7.0, 1.5) Å and exhibits a significantly large enhancement of γ/2 by a factor of ∼20 (γ/2 = 29.2 × 10 5 au) as compared to that in the S 0 state. The obtained γ value is found to be comparable to that in dicyclopenta[b;g]naphthaleno [1,2,3-cd;6,7,8-c′d′]diphenalene (NDPL) (γ = 22.8 × 10 5 au), 71 which is known to exhibit the largest possible TPA cross section (one of typical third-order NLO properties) in pure hydrocarbon systems of similar size. 22 In such molecular packing, the 1 (TT) state tends to be immediately created and to have a long lifetime, which indicates the possibility of utilizing the highly efficient thirdorder NLO property in the 1 (TT) state. The present results are expected to be applied to a variety of SF systems as a new design principle "SF-induced NLO" for the enhancement and the control of the NLO properties, which are predicted to be superior to those of conventional NLO systems. The investigations of a variety of SF-induced NLO systems are now in progress in our laboratory.

ACS Omega
Article ■ ASSOCIATED CONTENT

* S Supporting Information
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acsomega.9b02378.
Total Hamiltonian of the pentacene dimer model; molecular orbitals in the S 0 state; S 0 → S n and T 1 → T n excitation properties; accuracy of the excitation energies and NLO properties in the 1 (TT) state calculated by the UDFT method; numerical accuracy of the third-order NLO property using the FF method; molecular orbitals in the 1 (TT) state for A, B, and C;