Effects of Structural and Energetic Disorders on Charge Transports in Crystal and Amorphous Organic Layers

Understanding charge transports in organic films is important for both fundamental science and practical applications. Here, contributions of off-diagonal (structural) and diagonal (energetic) disorders to charge transports were clarified using molecular-based multiscale simulation. These disorders, important for understanding charge transport in organic systems, are investigated by comparing crystal and amorphous aggregates of N,N′-diphenyl-N,N′-bis(1-naphthyl)-1,1′-biphenyl-4,4′-diamine (NPD). Although NPD has been used as a hole transport material, it also exhibits comparable electron mobility experimentally. The experimental mobility and its electric field dependence in amorphous layers were reasonably reproduced by the multiscale simulation, confirming the electron transport properties of NPD. We assumed that the structural disorder would lower mobilities; however, the mobilities were found to be independent of the degree of structural disorder. Energetic disorder markedly lowered charge mobility instead. Charge migration in crystals was dominated by maximum electronic coupling pairs, whereas small electronic coupling pairs significantly contributed to charge transport in amorphous aggregate.

Understanding the charge transport properties of organic thin films is of great interest not only for fundamental science but also for practical applications in organic electronics such as organic light-emitting diodes (OLEDs), organic solar cells, and organic thin-film transistors 1 . Materials with high charge mobility have attracted much attention for use in these devices because they can greatly contribute to device performance. Organic molecules in devices are often in amorphous states, and charge transport in such systems is considered to occur by intermolecular charge hopping. Some models have been proposed to explain the charge transport in organic amorphous aggregates [2][3][4] . For example, the Gaussian disorder model 4 , one of the most well-known models for charge transport, describes the macroscopic charge transport process in amorphous organic solids. However, this model lacks a microscopic (molecular-level) description of charge transport processes. Recently, charge transport simulations based on Marcus theory 5 have been conducted for amorphous aggregates [6][7][8][9][10][11][12][13][14][15][16][17] , to reproduce experimental charge mobilities. However, the microscopic behaviour of charge transport in organic amorphous aggregates is still not fully understood.
Organic amorphous aggregates possess structural distributions by nature, such as distribution of intermolecular distance and that of relative orientation between adjacent molecules. The inherent structural distribution in an organic aggregate is usually called off-diagonal disorder or structural disorder. However, the effect of structural disorder on charge-hopping phenomena has not been fully considered at the molecular level. In this study, we perform charge transport simulations of crystalline and amorphous structures of N,N′-diphenyl-N,N′-bis (1-naphthyl)-1,1′-biphenyl-4,4′-diamine (NPD), a widely used hole transport material in OLEDs. A crystal and amorphous aggregate of NPD are used to exemplify ordered and disordered structures with negligible and large structural disorder, respectively. Energy levels of hopping sites are also distributed in amorphous aggregates; this distribution in energy is usually termed diagonal disorder or energetic disorder. We investigate the effect of energetic disorder as well as structural disorder on charge transfer through multiscale simulation by explicitly taking into account the organic molecule. Another purpose of this study is to examine the origin of bipolar charge

Results
For the NPD crystal, the structure determined by X-ray diffraction reported by Cheng et al. 32 was used without further structure optimisation (denoted as "cry-NPD", Fig. 1a). An amorphous aggregate of NPD was constructed by molecular dynamics (MD) simulation of a system composed of 4000 NPD molecules in a cubic cell (denoted as "amo-NPD", Fig. 1b). Figure 1b also exhibits three kinds of relative orientations in molecular pairs in amo-NPD as examples, showing various intermolecular packings. Supplementary Fig. S1 compares torsion angle distributions in amo-NPD and cry-NPD. They are significantly different due to the different nature of amorphous and crystalline aggregates. The reorganisation energies were calculated for molecules in the aggregated structure (λ agg + and λ agg − for hole and electron transport, respectively) by the quantum mechanics/molecular mechanics (QM/MM) method 33,34 . The electronic couplings for hole and electron transfer (H ij + and H ij − , respectively) were calculated according to ref. 35 . The computational details were also provided in Methods section. Using the values of H ij + , H ij − , λ agg + and λ agg − obtained above, the rate constants for hole and electron transfer (k ij + and k ij − , respectively) were calculated according to the following equation based on Marcus theory 5 , where T is the temperature, h is Planck's constant, k B is Boltzmann's constant, q is the charge of the carrier, F is the applied external electric field strength and Δx ij is the distance between molecules i and j along the direction of the applied electric field. ΔE i and ΔE j are the energies for the relevant two molecules i and j, respectively. On the basis of the method by Uratani et al. 17 , these energies were calculated with taking into account the difference in Coulombic interaction between neighbouring molecules. Figure 2 shows the highest occupied molecular orbital (HOMO) and lowest unoccupied molecular orbital (LUMO) of molecule I in cry-NPD (Fig. 1a). The HOMO and LUMO of molecule II were similar. Compared with the HOMO, which is relatively delocalised over the entire molecule, the LUMO is mainly localised on two outer naphthyl segments with only a small contribution along the central part, nitrogen atoms and biphenylene groups, as shown in Supplementary Table S1. A DFT-optimised NPD molecule used as the initial structure to construct amo-NPD provided a similar result (Supplementary Fig. S2 and Supplementary Table S1). The calculated H ij values for cry-NPD and amo-NPD are presented in Fig. 3. The maximum H ij values for electron (H ij − max = 8.8 meV for cry-NPD and 57.0 meV for amo-NPD) are larger than those for hole (H ij + max = 3.1 meV for cry-NPD and 11.9 meV for amo-NPD), suggesting that NPD possesses good electron transport properties. The calculated value of λ agg for electron transport is smaller than that of hole transport (Supplementary Table S2), which also indicates good electron transport properties. Taking a closer look at the case of cry-NPD reveals that the pair with the H ij − max of 8.8 meV has close face-to-face intermolecular contact between naphthyl segments ( Supplementary  Fig. S3a). In contrast, the molecular pair with H ij + max shown in Supplementary Fig. S3b has smaller intermolecular overlaps between orbitals (H ij + max = 3.1 meV). The H ij values of amo-NPD have a much larger distribution than those of cry-NPD (Fig. 3), and H ij − max and H ij + max of amo-NPD (57.0 and 11.9 meV, respectively; the molecular pairs are shown in Supplementary Fig. S4) are larger than those of cry-NPD. In the Gaussian disorder model 4 , H ij is solely determined by intermolecular distance. However, no clear correlation was found between H ij and intermolecular distance both for cry-NPD and for amo-NPD ( Supplementary Fig. S5). This result demonstrates that intermolecular orientational packing and frontier orbital distribution as well as intermolecular distance have significant influences on H ij .  Figure 4 shows the energetic disorders in cry-NPD and amo-NPD. The two crystallographically different molecules in the unit cell of cry-NPD, molecule I and II in Fig. 1a, have slightly different site energies from each other; the differences were 0.059 and 0.036 eV for electron and hole transports, respectively. Compared with that of cry-NPD, the energetic disorder in amo-NPD is much larger ( Fig. 4c and d), reflecting varied intermolecular interactions of surrounding molecules in the amorphous aggregate. In the Gaussian disorder model, the energetic disorder is assumed to be a Gaussian. However, in our work, no assumption was made for the distribution function. Figure 4c,d reveal that the distribution obtained in our calculation was slightly skewed to the higher  On the basis of k ij calculated by equation (1), we performed charge transport simulations by kinetic Monte Carlo method for both cry-NPD and amo-NPD. Figure 5a displays the dependence of the calculated hole (upper figure) and electron (lower figure) mobilities on the square root of electric field ( F) for cry-NPD with and without considering energetic disorder. The charge mobilities along the x-, y-and z-axes ( Supplementary Fig. S6 for respective axes) are different irrespective of the inclusion of energetic disorder, showing that charge transport is anisotropic in cry-NPD. The charge mobilities differ by one to two orders of magnitude depending on the direction of the applied electric field. Electron mobility is higher than hole mobility by one to two orders of magnitude, reflecting the larger H ij and smaller λ for electron transport than those for hole transport (see above). Reflecting the small energetic disorder in cry-NPD (Fig. 4a,b), slight decrease of charge mobility was found when energetic disorder is considered (Fig. 5a). Figure 5b shows the dependence of the hole and electron mobilities on F for amo-NPD. The difference in charge mobility along each axis is negligible for amo-NPD, reflecting the isotropic nature of the amorphous NPD system. Experimental charge mobilities are also presented in Fig. 5b. Our experimental data and those reported by Tse et al. 27 agree well with each other. For hole transport, positive F dependence, so-called Poole-Frenkel dependence, is experimentally observed. In contrast, electron mobility shows interesting F dependence; the mobility first decreases and then increases with increasing F. Compared with the experimental results, the charge mobility calculated without energetic disorder is one to two orders of magnitude higher at low F and has negative F dependence. When we consider energetic disorder, the hole mobility was one order of magnitude smaller than that obtained experimentally. In contrast, experimental electron mobilities are well reproduced by the multiscale simulation when the energetic disorder is considered. It should be emphasised that the experimentally observed F dependence was reproduced well for electron transport, including the negative F dependence at F lower than 300 V 1/2 cm −1/2 . From equation (1), F appears in the exponential term as λ In this case, k ij and therefore the velocity of electrons (v) do not strongly depend on F. As a result, charge mobility (which equals v/F) decreases with increasing F.
Comparison of Fig. 5a and b revealed that the calculated charge mobilities without considering energetic disorder and their F dependence in amo-NPD are comparable to those in cry-NPD (within the direction-dependent charge mobilities of cry-NPD). This holds for both hole and electron transfer, indicating that the structural disorder of NPD does not play a critical role in its charge mobility and applied electric field dependence. The result, "structural disorder in amorphous systems does not decrease charge mobility" is contrary to our expectation. Upon introduction of energetic disorder, a slight (one order of magnitude or smaller) decrease of charge mobility was found in the crystals, and the F dependence was mostly preserved. In contrast, for the amorphous NPD, the charge mobilities dramatically decreased when energetic disorder was considered, especially at low electric field, resulting in a positive F dependence, so-called Poole-Frenkel behaviour. These results show that it is not the structural disorder but the energetic disorder that causes the Poole-Frenkel behaviour of amorphous NPD. Hereafter, the molecular-level analysis will be carried out for electron transport because it reproduced the experimental results well. The same is essentially true for hole transport, although the simulated mobilities were about one order of magnitude smaller than experimental ones. Here, we introduce the migration distance for respective pairs, ΔX ij , as follows: ij ij diff where Δx ij is the distance between neighbouring relevant molecules i and j along the electric field. As shown later, charges hop in backward as well as forward directions. N diff is the number of difference of forward and backward charge hops between the molecules i and j. In this study, we calculated ΔX ij for the average of 10,000 trials of charge transport. Figure 6a,b show the correlation between ΔX ij and H ij − for cry-NPD with energetic disorder (the correlation for the case in which energetic disorder was ignored is summarised in Supplementary Fig. S7). For cry-NPD, large ΔX ij were obtained for the molecular pairs with maximum and second maximum H ij − . Figure 6c,d show the correlation between ΔX ij and H ij − for amo-NPD with energetic disorder. In this case, molecular pairs with small H ij − contribute markedly to charge transport (which provide large ΔX ij ), especially at F = 300 V 1/2 cm −1/2 , although negative values of ΔX ij are also found. The negative ΔX ij means that the molecular pair conveys charge in opposite direction to the electric field vector. The contribution of molecular pairs with small H ij − decreased at high electric field (Fig. 6d). Using ΔX ij , we can estimate why charge mobility in amorphous aggregates is comparable to that in crystalline aggregates when energetic disorder is not included. Figure 3 clearly shows that the H ij + and H ij − are widely distributed in amo-NPD. The distributions of H ij + and H ij − in cry-NPD are narrower and the maximum values are far smaller compared to the case of amo-NPD. In addition, Supplementary  Fig. S7c,d indicate that not only small H ij pairs but also large H ij pairs effectively convey charges in amo-NPD. To summarise, H ij values in amorphous aggregates are widely distributed due to the structural disorder, but almost the all distributed H ij pairs contribute the charge transport. Moreover, most of the H ij values in amo-NPD are larger than those in cry-NPD. This would be the reason that amorphous aggregates have good charge mobility comparable to crystalline systems under the condition that energetic disorder is not considered.
To understand the differences in charge hopping for cry-NPD and amo-NPD, we carried out further molecular-level analysis. Table 1 shows the number of charge hops in the forward direction (N fwd ) and backward direction (N bwd ), along with their sum (N all = N fwd + N bwd ) and difference (N diff = N fwd − N bwd ) during electron   100 nm. In the case of cry-NPD with energetic disorder, the pairs with the two highest H ij mostly contribute to charge hopping at 300 V 1/2 cm −1/2 (the sum of ΔX ij is 51.4 + 37.7 = 89.1 nm). For all the pairs, both forward (N fwd ) and backward (N bwd ) charge hoppings are found at 300 V 1/2 cm −1/2 ; the values of N all /N diff are 4.5-5.7, indicating that charge hops forward once every 4-6 hops. In contrast, at 1300 V 1/2 cm −1/2 , the N all /N diff values for the pairs with three highest ΔX ij values are almost 1, indicating that charges hop only in the forward direction at 1300 V 1/2 cm −1/2 . Table 1 also shows the case of amo-NPD. In amo-NPD, many molecular pairs are used during the charge transport and thus ΔX ij values for respective molecular pairs in amo-NPD are much smaller than those of cry-NPD. In amo-NPD, ΔX ij is relatively small for the pair with maximum electronic coupling (H ij − max = 57.0 meV) regardless of the inclusion of energetic disorder (smaller than 0.094 nm, Supplementary Figs S8a-S11a). The largest ΔX ij value (ΔX ij,max ) is determined by the pairs that have a smaller H ij compared to H ij − max (Supplementary Figs S8b-S11b). When the energetic disorder is considered, ΔX ij,max at 300 V 1/2 cm −1/2 is provided by the pair with H ij = 1.42 meV (Table 1 and Supplementary Fig. S8b). The second and third largest ΔX ij are provided by the pairs with H ij of 1.59 and 8.19 meV, respectively (Table 1 and Supplementary Fig. S8c). At 300 V 1/2 cm −1/2 , N fwd and N bwd are almost same order of magnitude for both cry-NPD and amo-NPD (Table 1). This result indicates that charge hopping not only in the forward and but also in the backward direction is common feature for both crystal and amorphous systems at low electric field. Finally, Fig. 7 compares electron transport trajectories in cry-NPD and amo-NPD considering energetic disorder (similar trajectories were found for hole transport). For the case of cry-NPD at 300 V 1/2 cm −1/2 (Fig. 7a), an electron reaches the counter electrode through mostly linear routes, with some small fluctuations; both forward and backward hops occur, as shown in Table 1. At 1300 V 1/2 cm −1/2 (Fig. 7b), an almost linear route is followed because of the repeated use of the molecular pairs with large H ij , although the route is tilted from the electric field direction (x-axis). In this case, backward hops seldom occurred, as mentioned above. In contrast, a much more complicated route is observed in amo-NPD at 300 V 1/2 cm −1/2 (Fig. 7c). Forward and backward charge hoppings occur in various directions, indicating diffusive-like behaviour in amo-NPD. At a higher applied electric field of 1300 V 1/2 cm −1/2 , the charge in amo-NPD reaches the counter electrode using more linear routes than that at 300 V 1/2 cm −1/2 (Fig. 7d).
In conclusion, we performed molecular-based multiscale charge transport simulations to investigate the effects of disorders in amorphous systems, in particular, the influence of structural disorder on charge transport in NPD, by comparing crystalline and amorphous aggregates. NPD has been used as a hole-transport material in OLEDs, but was found to exhibit electron as well as hole transport ability. Experimentally-obtained charge mobility and its field dependence in the amorphous thin layer, including interesting feature of the negative applied field dependence at low electric field, were quantitatively reproduced for electron transport by our multiscale simulation. The existence of "energetic" disorder markedly lowered charge mobility for amorphous NPD, especially at low electric field. However, the comparison of the crystalline and amorphous systems unexpectedly revealed that "structural" disorder does not obviously decrease charge mobility. Analysis of the actual contribution of each molecular pair to the migration distance revealed a substantial difference between crystalline and amorphous systems. Charge migration was basically dominated by the molecular pairs with maximum and second maximum electronic coupling in the crystalline aggregate, whereas the molecular pairs with relatively small electronic coupling contributed to the charge transport in the amorphous aggregate. These findings on the effects of structural and energetic disorders on charge-hopping phenomena at the molecular level should aid the design of materials with superior charge transport properties.

Methods
Computational details. The MD simulations of amorphous systems were performed using the LAMMPS program 36 . For the MD calculation, the Dreiding force field 37 was used and the bond lengths of the DFToptimised molecular structure were applied as the equilibrium bond length parameters. The DFT calculation (B3LYP functional 38 and 6-31G(d) basis set) was performed on the Gaussian 09 program package 39 . To mimic the deposition process, the initial MD simulation was performed in the constant-volume, constant-temperature (NVT) ensemble at a density of 0.1 g cm −3 and temperature of 573 K for 10 ps. Then, the MD simulation was performed under the constant-pressure, constant-temperature (NPT) ensemble for 1.0 ns at 298 K until the density of the system became constant. Finally, the structure was geometry optimised to eliminate the structure deviation originating from vibrational motion during MD calculations. The final density of the amorphous structure was 0.95 g cm −3 (the density of cry-NPD was 1.22 g cm −3 ). The amorphous structure thus obtained, amo-NPD, is shown in Fig. 1b. λ agg + and λ agg − were calculated by QM/MM method 33,34 , which includes the effect of intermolecular interactions 17, 40 . In the QM region, DFT with B3LYP functional and 6-31G(d) basis set was used and the Dreiding force field were used in the MM region for surrounding molecules within a distance of 30 Å. For comparison, the reorganisation energies for an isolated NPD molecule (λ iso ) were also calculated according to previous reports 16,17,41,42 . The values of λ agg and λ iso are summarised in Supplementary Table S2.
H ij + and H ij − were calculated by the following equation.
Here, Ĥ ij is the electronic Hamiltonian of a dimer system composed of molecules i and j. ψ i and ψ j are the HOMO or LUMO of the isolated molecules i and j, respectively. Molecular pairs with a centre-to-centre distance within 25 Å are considered in the H ij calculation.
The charge transport simulations for both cry-NPD and amo-NPD were performed based on kinetic Monte Carlo calculations. The distance (L) of charge migration along the applied electric field F is 100 nm and 10,000 trials of charge transports were simulated to calculate the travelling time (t). The mobility μ was calculated from μ = L/(tF), where F = |F|. Data availability. The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.