A general charge transport picture for organic semiconductors with nonlocal electron-phonon couplings

The nonlocal electron-phonon couplings in organic semiconductors responsible for the fluctuation of intermolecular transfer integrals has been the center of interest recently. Several irreconcilable scenarios coexist for the description of the nonlocal electron-phonon coupling, such as phonon-assisted transport, transient localization, and band-like transport. Through a nearly exact numerical study for the carrier mobility of the Holstein-Peierls model using the matrix product states approach, we locate the phonon-assisted transport, transient localization and band-like regimes as a function of the transfer integral (V) and the nonlocal electron-phonon couplings (ΔV), and their distinct transport behaviors are analyzed by carrier mobility, mean free path, optical conductivity and one-particle spectral function. We also identify an “intermediate regime” where none of the established pictures applies, and the generally perceived hopping regime is found to be at a very limited end in the proposed regime paradigm.


I. THE MATRIX PRODUCT STATES PARAMETERS
In the following we present details on our MPS time evolution parameters and respective numerical convergence verification. The calculation of C(t) involves three sets of time evolution: the imaginary time evolution for Ψ β , real time evolution for e −iĤt |Ψ β and real time evolution for e −iĤtĵ (0) |Ψ β . For efficient calculation, we use different bond dimensions for the three time evolutions. More specifically, time evolutions for Ψ β and e −iĤt |Ψ β share the same bond dimension iM, and the time evolution for e −iĤtĵ (0) |Ψ β has the bond dimension M. Other important parameters relevant to the time evolution include the number of molecules in the periodic one dimensional chain N, the number of truncated harmonic oscillator eigenbasis for the intermolecular vibration l, and time evolution step dt. For most of our calculations, iM, M, N, l and dt are set to 64, 80, 21, 40, 50 respectively, and numerical convergence in the large transfer integral V = 150 meV and strong nonlocal EPC strength ∆V = 60 meV limit is shown in Fig. 1. This is considered to be a set of challenging parameter due to the strong EPC and non-trivial physical behavior as discussed in the main text. We note, however, that a fixed set of MPS parameters is not sufficient for converged result over the whole (V, ∆V ) plane explored in the main text, and we in fact adjust the MPS parameters accordingly to ensure that nearly exact result is obtained. We show in Fig

II. BENCHMARK OF THE ONE-PARTICLE SPECTRAL FUNCTION
The evaluation of the one-particle spectral function A(k, ω) further takes advantage of the fact that the thermal equilibrium state is a zero-electron state. In this case, the thermal field dynamics algorithm can be reformulated to reduce computational cost 1,2 . For the calculation of A(k, ω) we

III. THE NUMBER OF INTRAMOLECULAR VIBRATION MODES
Organic molecules usually contain dozens of atoms and it is common for more than 10 vibrational modes to contribute to local EPC. For simplicity, theoretical studies often reduce the modes into one effective mode 3,5,6 . However, by taking this approach, the correlation function C(t) typically exhibits spurious long correlation time which is generally not possible for realistic material. In such cases, an artificial and sometimes arbitrary broadening is usually applied to C(t) in order to mimic realistic world and determine an absolute value for mobility. In the main text most of our results are reported with a 4-(intramolecular-)mode model and C(t) generally rapidly decays to zero. But in cases of small V and ∆V we have to further resort to a 9-(intramolecular-)mode model with vibration energy and EPC constant listed in Table I to avoid artificial recurrence.
Starting from more than 40 vibrational modes that contribute to local EPC in rubrene crystal obtained by ab initio DFT calculation, the 9-mode model is generated by dropping the modes with λ m = g 2 m,I ω m < 20 cm −1 and adding λ m of the dropped modes to the closest retained modes, while the 4-mode model is generated in the same way except that the dropping threshold is upraised to 50 cm −17 . Thus, the 9-mode model and the 4-mode model share the same total reorganization energy ∑ m λ m and the 9-mode model in principle describes realistic material better than the 4 mode-model.
In Fig. 5 we list the MPS parameters for which we have used the 9-mode model as well as the corresponding correlation functions C(t). In Fig. 5(a) and (c) it can be seen that that the artificial recurrence around 6000 a.u. for the 4-mode model is absent in the 9-mode model.