Predicting Phosphorescence Quantum Yield for Pt(II)-Based OLED

Subscriber access provided by CLARKSON UNIV

C: Plasmonics; Optical, Magnetic, and Hybrid Materials

Predicting Phosphorescence Quantum Yield for Pt(II)Based OLED Emitters from Correlation Function Approach Inkoo Kim, Won-joon Son, Youn-Suk Choi, Alexey A. Osipov, Dongseon Lee, Hasup Lee, Yongsik Jung, Jhunmo Son, Hyeonho Choi, Wataru Sotoyama, Alexander Yakubovich, Jaikwang Shin, and Hyo Sug Lee J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.9b02031 • Publication Date (Web): 05 Apr 2019 Downloaded from http://pubs.acs.org on April 5, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Predicting Phosphorescence Quantum Yield for Pt(II)-Based OLED Emitters From Correlation Function Approach Inkoo Kim,† Won-Joon Son,∗,† Youn-Suk Choi,∗,† Alexey Osipov,‡ Dongseon Lee,† Hasup Lee,† Yongsik Jung,† Jhunmo Son,† Hyeonho Choi,† Wataru Sotoyama,¶ Alexander Yakubovich,‡ Jaikwang Shin,† and Hyo Sug Lee† †Samsung Advanced Institute of Technology (SAIT), Samsung Electronics, 130 Samsung-ro, Yeongtong-gu, Suwon 16678, Korea. ‡Samsung R&D Institute Russia (SRR), Samsung Electronics, 12 Dvintsev Street, Moscow 127018, Russia. ¶Samsung R&D Institute Japan (SRJ), Samsung Electronics, 2-7 Sugasawa-cho, Tsurumi-ku, Yokohama 230-0027, Japan. E-mail: [email protected]; [email protected] Abstract We report a new formulation for Golden Rule based predictions of photoluminescence quantum yields (PLQY) of phosphorescent emitters containing a heavy element, and its implementation compatible with first-principle computation frameworks. The main components of PLQY (i.e., phosphorescent rate and inter-system crossing rate from the lowest triplet state to the ground singlet state) are obtained through correlation functions in time domain, and the relativistic effects are also considered using the relativistic effective core potentials. The spin–orbit coupling is treated in a perturbative manner to generate spin–orbit corrected, two-component T1 substates within

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 35

single-excitation theory, where the electronic transition dipole moments and the nonBorn-Oppenheimer coupling matrix elements to the S0 state are computed. We applied this new approach to photophysical properties of 34 Pt(II) complexes designed for the organic light-emitting diode (OLED) applications, and observed a good agreement between predictions and experiments over diverse scaffolds. For the two representative complexes, further analysis on the nonradiative characteristics was performed based on the decomposition of the non-Born–Oppenheimer coupling into contributions from the nuclear vibrations and from the excited-state electronic structure.

1

Introduction

Photoluminescence quantum yield (PLQY), characterizing the radiative ability of light emitters, 1 is a representative photophysical property often expressed as

ΦPL =

kr , kr + knr

(1)

with the radiative decay rate kr and the nonradiative relaxation rate knr . Presuming Kasha’s rule, two types of photoluminescence can be considered: fluorescence of radiative transition from excited singlet states manifold, and phosphorescensce of radiative transition from excited triplet states manifold. The most well-recognized and, arguably, most important industrial application of phosphorescence can be found in organic light-emitting diodes (OLED), 2 encountered in high-value display applications ranging from smartphones to TV panels. 3,4 High-performance OLED devices based on all-phosphorescent emitters have been de facto the most promising candidates for next-generation displays, due to their high contrast ratios, wide color gamut and power efficency. 5 During the last decade, the importance of cyclometallated Pt(II) and Ir(III) complexes as light emitters has become firmly entrenched, because of their photo- and electrochemical stability, tunable emission color and especially their high PLQY. 6

2


Page 3 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


In OLED device, external electroluminescence quantum efficiency (EQE) is normally regarded as the most important device parameter, where the PLQY of the emitter is inherently related. 7 The EQE, in short, describes a ratio between the number of emitted photons and the injected charge carriers, expressed as the product of light extraction efficiency ηout , upper limit fraction of radiatively decaying excitons ηS/T , charge-carrier balance factor γ, and the PLQY, i.e., ηEQE ≡ ηout · ηS/T · γ · ΦPL . For emitters containing a heavy transition element, ultrafast ISC between excited singlet states and triplet states (especially, the T1 state) occurs, leading the ηS/T factor close to unity. 8,9 Sophisticated device optimizations can also contribute to enhancing EQE, for instance, by optimizing ηout and γ, 10 but a deliberate design of emitters with high PLQY is generally considered as on of the most critical issues in the development of novel OLED systems. Theoretical prediction of the phosphorescence quantum yield, however, has been a notoriously difficult task, since comprehensive theoretical description of the underlying process of radiative and nonradiative triplet-to-singlet transitions demand the intertwined and simultaneous treatment of the spin–orbit coupling and the nonadiabatic effects. 9,11 The theoretical efforts for the first-principles approaches have still been continued, and recent reports have made remarkable advances. Among them, one can argue that time-dependent approach based on the correlation function formalism 12,13 in the framework of Golden Rule expression 14 appears to be promising for the rate calculations as it provides an analytic path-integral formulation under the harmonic oscillator approximation. Using this approach, Marian and co-workers 15,16 derived the rate expressions for the inter-system crossing process based on the first-order vibronic spin–orbit coupling. Shuai and co-workers 17 further extended the approach to include the spin-vibronic perturbation terms arsing from the second-order perturbation treatment of the spin–orbit and the nonadiabatic vibronic couplings. This approach has also been applied to simulate the optical spectra and fluorescencent rates of various molecules, 18–20 and recently, for phosphorescent rates of purely organic molecules. 21 For a rigorous description of ISC, one should construct wavepacket propagation on spin–orbit cou-

3



pled potential energy curve, 22,23 however, this is generally not feasible for large molecules due to enormous computational cost. Escudero 24 reported a computational methodology for predicting PLQYs by constructing the deactivation kinetics of Ir(III)-based emitters. 25 By incorporation of kinetic model into the time-dependent approach, Zhang et al. 26 recently simulated the temperature dependence in PLQY for three Ir(III) complexes. The applicability and validity of new methods are, however, often exemplified by series of analogous structures with a restricted set of functional groups. Timely development of a method for predicting PLQY reliable over wide range of materials is still in need, especially for effective computational screening of novel phosphorescent emitters in industries. In this work, we have devised a quantitative method enabling such predictions based on the correlation function approach, 13,18 with which the Golden Rule equations for the phosphorescence rate and the inter-system crossing rate are solved analytically. from the excited T1 to the ground S0 state. The spin–orbit coupling is efficiently treated in the first-order perturbation theory using the relativistic effective core potential (RECP). 27 Applications are made to the photophysical properties for the Pt(II) complexes designed for the OLED applications, 28 and good agreement with experiments was observed.

2

Theory

All equations are written in atomic unit.

2.1

Transition rates

The electronic T1 state of total spin of S = 1 is characterized by three degenerate substates (M = 0, ±1). At room temperature, fast thermalization of the substates allows the transition rate (either radiative or nonradiative) originating from the T1 state to be approximated as

4


Page 4 of 35

Page 5 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


an average of the substates:

k=

1X M k . 3 M

(2)

The energy splittings between these substates for usual phosphorescent emitters incorporating a 5d-block transition element are order of a few hundred wavenumbers, and hence sufficient thermal energy for equally populating these substates should be available at room temperature. 29 In what follows, the so-called high temperature limit is assumed; in the case of low temperature, one should consider the thermal distribution of the substates.

2.2

Emissive triplet state and spin–orbit correction

In the single-excitation theories, such as linear-response time-dependent density functional theory (TDDFT) 30,31 within Tamm-Dancoff approximation 32 or configuration interaction singles (CIS), the nonrelativistic T1 substates can be expressed as 11 1 X |T10 i = √ Cia (|φia i + |φ¯i¯a i) , 2 i,a √ X Cia |φ¯ia i , |T11 i = 2

(3) (4)

i,a

|T1−1 i =

√ X 2 Cia |φi¯a i ,

(5)

i,a

where Cia is the transition density matrix element corresponding to a singly excited configuration |φia i from i-th occupied to a-th virtual orbital. Here, the unbarred and barred indices represent α- and β-spin orbitals, respectively. In practice, only the |T10 i wavefunction is obtained as the solution of eigenvalue problem, and the other two substates (M = ±1) can be generated by applying the spin-ladder operators on the |T10 i. These states define the orthonormal, spin-adapted zeroth-order states subject to spin–orbit correction for describing the emissive T1 state of the phosphorescent emitters.

5



Page 6 of 35

As a result of strong spin–orbit coupling of the central heavy element, the T1 substates are significantly perturbed by nearby singlet excited states, thereby the spin-forbidden transition to the S0 state can be allowed. 33 In the first-order perturbation theory, the two-component, ˆ SO can be written as 34,35 perturbed T1 substate by the spin–orbit Hamiltonian H |T˜1M i = |T1M i +

X

|Sn i

n=1

ˆ SO |T M i hSn |H 1 , ET1 − ESn

(6)

where the tilde represents the spin–orbit perturbed state, ET1 and ESn are the energies of the T1 and the n-th excited singlet Sn state, 1 X Sn |Sn i = √ Cia (|φia i − |φ¯i¯a i) , 2 i,a

(7)

respectively. Similar to eq 6, the S0 state may also be written in a perturbative manner with the triplet manifold; however, by noting that the T1 –S0 gap is at least 2 eV and above for visible photons, one may neglect triplet contributions, 36 hence the following approximation can be adopted:

|S˜0 i ≈ |S0 i .

(8)

The orthogonality condition hS0 |T˜1M i = 0 is satisfied for all M , due to orthogonality of singlet and triplet manifolds. The non-vanishing electronic transition dipole moment for phosphorescence can then be written as

hS0 |ˆ µ|T˜1M i =

X hS0 |ˆ ˆ SO |T M i µ|Sn i hSn |H 1

ET1 − ESn

n=1

,

(9)

where µ ≡ −eˆ r is the electronic dipole operator. This recapitulates the well-known intensity borrowing effect from the excited singlet states. 8 One should note that the approximation eq 8 eliminates the need for computing the transition dipole moment between excited triplet

6


Page 7 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


states, which would require, for instance, in TDDFT the quadratic-response formalism. 37,38

2.3

Spin–orbit couplings

We are now to compute the spin–orbit coupling between the electronic states, entering in eq 6. The spin–orbit coupling of heavy elements is most conveniently described within the two-component framework of spin–orbit RECP, containing the (spin–orbit-)averaged RECP (AREP), which replaces the core electrons and also incorporates the scalar relativistic effects (mostly from mass-velocity and Darwin terms) and the effective one-electron spin–orbit operator, and is given by 27

ˆ SO , Uˆ RECP = Uˆ AREP + H

(10)

with

Uˆ AREP = ULAREP (r) +

l L−1 X X l=0 m=−l

× [UlAREP (r) − ULAREP (r)] |lmi hlm| , ˆ SO = sˆ · H

L X l=1

×

2 ∆UlRECP (r) 2l + 1

l l X X m=−l

(11)

|lmi hlm|ˆl|lm0 i hlm0 | ,

(12)

m0 =−l

RECP RECP RECP where UlAREP (r) = (2l + 1)−1 [l · Ul,l−1/2 (r) + (l + 1) · Ul,l+1/2 ] and ∆UlRECP (r) = Ul,l+1/2 (r) − RECP Ul,l−1/2 (r). Here, l and m denote the orbital angular momentum quantum number, and the

magnetic quantum numbers for the given l, and |lmi hlm| is the two-component projection operator. Note that the first summations in eqs 11 and 12 are constrained to a certain angular momentum L at which the contribution from the higher angular momenta are summed into a single term or neglected entirely. The Uˆ AREP term enters in the definition of the nonrelativistic electronic Hamiltonian and 7



Page 8 of 35

ˆ SO are transformed to can be treated self-consistently. The atomic spin–orbit integrals of H molecular orbital basis hSO , and then used for calculating the spin–orbit matrix elements ˆ SO |T1M i. Note that, between the singly-excited configurations that are required for hSn |H first, since the one-electron spin–orbit operator is employed, the matrix elements vanishes when two or more indices differ between the configurations, and secondly, the transition density matrix of the Sn state eq 7 has a spin-block diagonal structure, and therefore, only the appearance of the hφia | and the hφ¯i¯a | configurations are expected. Then, for the given i-th closed orbital,

ˆ SO |φib i = hSO , hφia |H ab

(13)

ˆ SO |φ¯ib i = hSO hφia |H a¯b ,

(14)

ˆ SO |φ¯ib i = hSO hφ¯i¯a |H a ¯b ,

(15)

ˆ SO |φ¯i¯b i = hSO hφ¯i¯a |H a ¯¯b ,

(16)

and for the given a-th virtual orbital,

ˆ SO |φja i = −hSO , hφia |H ji

(17)

ˆ SO |φ¯ja i = −h¯SO , hφia |H ji

(18)

ˆ SO |φj¯a i = −hSO hφ¯i¯a |H j¯i ,

(19)

ˆ SO |φ¯j¯a i = −h¯SO hφ¯i¯a |H j¯i ,

(20)

are the eight-types of the surviving integrals. The time-reversal symmetry constraints, which replaces the nonrelativistic spin symmetry in the relativistic regime, leads to the following relations with general orbital indices p and q: 39,40

SO∗ hSO pq = hp¯q¯ ,

(21)

SO∗ hSO p¯q = −hp¯ q .

(22)

8


Page 9 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


With this in hand, the spin–orbit coupling in eq 6 can now be evaluated by using eqs 13–20, and we finally arrive at the following expressions:

ˆ SO |T 0 i = i hSn |H 1

nX

Sn Cia

X

i,a

+ (−1)

CibT1 Im[hSO ab ]

b

i+j+1

X

Sn Cia

X

i,a

ˆ SO |T1−1 i = hSn |H

X

Sn Cia

X

i,a

T1 Cja

Im[hSO ji ]

o ,

(23)

j

CibT1 hSO a¯b

b

+ (−1)i+j+1

X

Sn Cia

i,a

X

T1 SO Cja hj¯i ,

(24)

j

ˆ SO |T11 i = hSn |H ˆ SO |T1−1 i∗ . hSn |H

(25)

ˆ SO has led Note that, in the derivation, the time-reversal symmetry in the construction of H to eq 25.

2.4

Phosphorescence rate

The nonvanishing transition dipole moment in eq 9 due to the spin–orbit coupling leads to the spontaneous transition from the excited T1 state to the electronic ground S0 state accompanied by a photoemission. Under the dipole approximation with the minimal coupling Hamiltonian, 41 the high-temperature T1 →S0 phosphorescence decay rate can be approximated as the average over substates with different spin-projections M of the integrals of the

9



Page 10 of 35

differential photon emission rate, i.e., Z 2 1 X ∞ 4ω 3 kphos = dω 3 hS0 |ˆ µ|T˜1M i 3 M 0 3c X 2 × Pν (T ) |hν 0 |νi| δ(∆E + Eν − Eν 0 − ~ω) ν,ν 0

2 Z ∞ X 4ω 3 2 ˜ dω Pν (T ) |hν 0 |νi| µ|T1 i = 3 hS0 |ˆ 3c 0 ν,ν 0 × δ(∆E + Eν − Eν 0 − ~ω),

(26)

where ν, ν 0 label the vibrational states of the T1 and the S0 energy levels, respectively, with Eν , Eν 0 being their respective vibrational energies, ∆E is the adiabatic energy difference, and Pν = exp(−Eν /kB T )/Z describes the Boltzmann distribution for the initial T1 vibrational manifolds at temperature T where Z being the vibrational partition function.

The substate-averaging can be effectively included within the prefactor (i.e., P hS0 |ˆ µ|T˜1 i = 13 M hS0 |ˆ µ|T˜1M i). Note that with eq. 25, we can avoid the explicit calculation of hS0 |ˆ µ|T˜1−1 i. The appearance of the delta function represents the energy conservation law. We have assumed that the vibrational manifolds of all T1 substates are identical within the high-temperature approximation. The computation of eq 26 involves a large number of calculations of the Franck-Condon factors hν 0 |νi from all combinations of the vibrational states of the initial T1 and final S0 states. Thus, one should resolve the problem of an efficient calculation of the integrand of eq 26 for practical applications. This is dealt in section 2.6 where the correlation function approach is introduced for the analytical solution of this problem.

10


Page 11 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2.5

Inter-system crossing rate

In the high-temperature limit, nonradiative T1

S0 ISC rate can be determined using Golden

Rule expression 14

kisc =

2 1 X 2π X ˆ nBO |T˜M , νi Pν (T ) hS0 , ν 0 |H 1 3 M ~ ν,ν 0

× δ(∆E + Eν − Eν 0 ),

(27)

ˆ nBO describes the breakdown of the Born– where the non-Born–Oppenheimer Hamiltonian H Oppenheimer adiabatic approximation. We note that the interaction Hamiltonian used here ˆ SO was also included in the interaction Hamiltonian. Our apdiffers from ref 17, where H proach can be regarded as a two-step approach for the calculation of the nonadiabatic coupling between different spin states. The approximate two-component relativistic wavefunctions by the first-order (spin–orbit-)correction are generated before the calculation of the nonadiabatic coupling. This two-step approach should work whenever the states of interest, namely the S0 and the T1 , are well separated in vicinity of the equilibrium geometries so that eqs 6 and 8 are satisfied. The coupling term in eq 27 can be expressed in the perturbative manner similarly to eq 9 but, this time, with vibronic states as the nonadiabatic transition is mediated by the momentum operator acting on the nuclear coordinates: 2 ˆ nBO |T˜M , νi hS0 , ν 0 |H 1 2 X X hS , ν 0 |H ˆ nBO |Sn , ν 00 i hν 00 |νi hSn |H ˆ SO |T1M i 0 = , n 00 ET1 − ESn

(28)

ν

where ν 00 and ESn denote the vibrational state and the energy of the Sn state, respectively. The vibronic state can be decomposed into the direct product of the electronic and the vibrational states, for instance, |T1M , νi = |T1M i ⊗ |νi. We note that the nuclear coordinate

11



Page 12 of 35

ˆ SO has been neglected as well as the vibrational energies in the denominator, dependence of H the latter of which is known as Plazcek approximation. Under harmonic and Condon approximations, the non-Born–Oppenheimer coupling matrix elements between states of the same spin symmetry can be estimated as 42

ˆ BO |Sn , ν 00 i ≈ −~2 hS0 , ν 0 |H

X

hSn |

κ

∂S0 ∂ν 0 i hν 00 | i, ∂Qκ ∂Qκ

(29)

where Qκ is the κ-th normal coordinate of either the T1 or the S0 state. The normal coordinate systems of different electronic states can be used interchangeably as each set of coordinates describes the complete nuclear motion near their equilibrium geometries. SubP stituting eq 29, and using the identity ν 00 |ν 00 i hν 00 | = 1, and defining the vector

tM κ

= −i~

∂S0 ˆ SO |T M i X hSn | ∂Q i hSn |H 1 κ

ET1 − ESn

n

,

(30)

we rewrite eq 28 as 2 X 0 ∂ν 0 0 ˆ BO ˜ M M M ∗ ∂ν i. tk tl h |νi hν| hS0 , ν |H |T1 , νi = ∂Qκ Qκ0 κ,κ0

(31)

Additionally, the averaging over the substates in eq 27 can be carried out by defining the matrix

T= =

1 X M M† t ·t 3 M 1 Re[T0 + 2T1 ], 3

(32)

which is real and symmetric by construction. This can be easily shown by examining ˆ SO |T10 i is purely imaginary, and hSn |H ˆ SO |T11 i and eqs 23–25; the spin–orbit integral hSn |H ˆ SO |T1−1 i are related through time-reversal symmetry. Using eqs 31 and 32, kisc is finally hSn |H

12


Page 13 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


expressed as

kisc =

X 2π X ∂ν 0 ∂ν 0 Tκκ0 Pν h |νi hν| i ~ κ,κ0 ∂Q Q κ κ0 0 ν,ν

× δ(∆E + Eν − Eν 0 ).

(33)

The form of eq 33 is identical with that derived for the IC process 13 in which Tκ,κ0 will be replaced by the nonadiabatic coupling matrix elements between states of the same spin. This should be expected since, in the relativistic regime when the spin is no longer a good quantum number, the ISC and the IC processes cannot be differentiated from each other.

2.6

Correlation function approach

The Fourier transformation of the delta function, i.e., δ(E) = (2π)−1

R∞ −∞

dτ exp(iEτ ) with

τ ≡ t/~, to the rate expressions of eqs 26 and 33 yields 12,13,17

kphos

Z ∞ 2 Z ∞ 2ω 3 ˜ = µ|T1 i dω dτ hS0 |ˆ 3πc3 Z 0 −∞ × exp[−i(~ω − ∆E)τ ]

0 ˆ HO ˆ HO )], × Tr[exp(−iτ 0 H ) exp(−iτ H Z ∞ 1 X 0 kisc = Tκκ dτ ~Z κ,κ0 −∞

ˆ 0 )Pˆκ0 exp(−iτ H ˆ HO )], × Tr[Pˆκ exp(−iτ 0 H HO

(34)

(35)

ˆ HO and H ˆ 0 are the multiple harmonic oscillator Hamiltonian where τ 0 ≡ −τ −i/(kB T ), and H HO ˆ HO |νi = Eν |νi with Eν = (k + 1 )~ωk for the vibrational states of the T1 and S0 states that H 2 where k is the vibrational quantum number, and Pˆκ = −i~∂/∂Qκ is the momentum operator in the normal coordinates. We reiterate that the substate-averaging has been accounted for by averaging the necessary coupling factors in hS0 |ˆ µ|T˜1 i and also in Tκκ0 . The correlation functions given by the trace functions above can be solved analytically in the path integral 13



Page 14 of 35

formalism at the given moment of time t, thus rendering direct solutions of the kphos and the kisc via the back Fourier transformation. The algorithm details in the unified time-integration scheme for these correlation functions are succinctly described in refs 13, 18 and 20. The correlation function approach described above has been implemented in a selfcontained Fortran code. The code can also perform a normal mode analysis including Duschinsky transformation in terms of curvilinear internal coordinates, 43 which is crucial since rectilinear internal coordinate often overestimate the displacement and consequently, the reorganization energy. In addition, we employed the displaced-distorted independent harmonic oscillator approximation 44 in which the potential energy surfaces (PES) for the states of interest are described by two sets of harmonic surfaces, with shifted equilibrium positions and different frequencies, and thereby excluding the mode-mixing effect (which corresponds to rotation of the one PES relatively to the other) from our model. The normal modes of the S0 state modes are rearranged such that a maximum overlap with the modes of the T1 state is achieved, and then the Duschinsky matrix is set to unity:

QS0 = QT1 + d,

(36)

where QS0 and QT1 are the respective sets of normal coordinates of the S0 and the T1 states, and the vector d represents the displacement between the two minima of potential energy surfaces. For large molecules such as phosphorescent emitters, a significant change in molecular geometries can lead to inaccurate representation of curvilinear Duschinsky matrix, and its determinant often deviates from unity. 15 Such deviations are undesirable and should be eliminated prior to the integration of correlation functions in time-domain. Therefore, we have taken the simplest approach for this problem by undertaking the displaced and distorted harmonic oscillator model, which will be shown later that the difference between the rates with and without the mode-mixing effect are marginal. Next, we move on to discuss the replacement of the delta function in the rate expres-

14


Page 15 of 35

16526

b)

a)

9

kisc (s−1)

10 Time (fs)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


40 20 0 -3

3 0 10 ) 0 0 1 ( s 3 -3 Real a axi xis (10 10 ag ) Im

τL=100 fs

108 107 6

1.89 x 106 s−1

5

2.08 x 105 s−1 3.93 x 104 s−1

10 10

104 0

10000

20000

Energy (cm−1)

1000 fs

10000 fs

Figure 1: (a) Correlation functions for kisc with projections of real and imaginary parts for various τL , and (b) the corresponding Fourier-transforms in frequency domain. The adiabatic energy difference is drawn as a dashed line. sions by the Lorentzian line shape function, L(E; Γ) =

Γ 1 π E 2 +Γ2

in frequency-domain or

L(t; τL ) = exp(−τL−1 t) in time-domain, with full-width at half-maximum (FWHM) of 2τL −1 . Such substitution does not only ensure the convergence of the time integration, but when a suitable damping parameter τL is chosen, it also partially compensates the limitations of the harmonic approximation by energy level broadening. An example of the time evolution of the ISC correlation function in eq 35 is shown in Figure 1a. The effect of the damping parameters between 100 fs through 10000 fs are almost indistinguishable in time-domain, however, these parameters determine the asymptotic behavior of the functions in frequency domain (Figure 1b) after the back-transformation. The selection of decay parameter is a subtle issue and in this work we obtain this value by comparing the results with respect to the experiments. We note that, as can be seen in Figure 1b, the kisc is maximized near the reorganization energy and then decays exponentially, demonstrating the energy gap law. The most time-consuming part is the construction of the correlation functions in eqs 34 and 35. We have mitigated the demanding computations by a simple and straightforward parallelization of the code, thanks to the fact that each time step is independed from the previous one. The time integration interval of [−6553.6:6553.6] fs with a time step of 0.1 fs, was found to be satisfactory on both cost and precision, and hence was used for all correlation functions. The Fourier transformation in one-dimension was computed using FFTW library. 45 Lastly, the computational cost of the correlation function approach is negligible compared 15



to that of Hessian or excited states calculation. The largest molecule considered in this work is the complex 21 (introduced in section 3) containing 113 atoms, and the corresponding kisc calculation took about 500 seconds using 32 CPU cores of two Xeon E5-2698 v3 at 2.30 GHz, while the numerical calculations of normal modes with 678 single calculations of displaced geometries in the same compute node would be roughly 60 times longer.

3

Applications

Cyclometallated Pt(II) complex with a tetradentate ligand has recently emerged as an important class of phosphorescent emitters for OLED devices. 46 Tetradentate ligands, often considered as good chelators, can provide rigid scaffolds, rendering stiffness of molecular structure upon excitations. Stable Pt(II), a formally d8 ion, is known to prefer a squareplanar structure with an empty dx2 −y2 orbital of a very strong anti-bonding character. 47 Stability and photoluminescence efficiency of the complex is closely related to the population of this formally empty orbital, induced either by photo- or electro-excitation, or even by thermal population. The anti-bonding character of dx2 −y2 weakens metal-ligand bonds and introduces a significant distortion of the square-planar geometry. These factors unlock an additional efficient nonradiative decay channel, deteriorating the PLQY of emitters. 48 In this respect, the rigid scaffold provided by tetradentate ligands should prevent the distortion of the Pt-ligand bonds, and hampers the population of the metal-centered states. Their radiative ability of the T1 state is, therefore, eventually governed by the interplay of spin–orbit and nonadiabatic couplings to the ground S0 state as described in eqs 26 and 27. We applied the present method for predicting the PLQY of 34 Pt(II)-based phosphorescent OLED emitters in the visible color range (for chemical structures, see Supplementary Information Figure S1). Experimental data have been adopted from ref 28, and the original numbering scheme for the complexes has been retained here. The set contains various tetradentate coordination modes, thus, ensuring the diversity of both ligand structures and

16


Page 16 of 35

Page 17 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


PLQY values in order to minimize bias in assessing the performance of our model.

3.1

Computational details

All DFT calculations were performed with the Turbomole program suite. 49 The molecular geometries were optimized by DFT using the hybrid version of generalized gradient approximation (GGA) functional of Perdew–Burke–Ernzerhof 50,51 with the resolution-of-the-identity approximation. 52 All-electron def2-SVP basis sets 53 were used for all the atoms, except the Pt atom, which was treated using the triple-ζ quality basis sets and a 60-electron core quasi-relativistic effective core potential. 54 The spin-unrestricted formalism was employed in geometry optimization of the T1 state to calculate the adiabatic energy difference ∆E to the S0 state. A single-point restricted DFT calculation was also performed at this geometry to generate the reference function for time-dependent DFT calculations under Tamm-Dancoff approximation 32 for vertical excitation energies of the 50 excited singlet and the T1 state as well as nonadiabatic coupling matrix elements in Cartesian coordinates, 55 which are then transformed to normal coordinates using the results of the harmonic force constant calculations. Solvent effects on the ground-state geometries and excitation energies were included using the conductor-like screening model (COSMO) of the dielectric continuum solvation method with default parameters. 56 Numerical frequency calculations were carried out to obtain the normal modes and ensure the minimum energy geometry. The optimized geometries and the corresponding vibrational frequencies are provided in Supplementary Information Tables S1 and S2. We note that the complexes 24, 36–39 are omitted from the original set as the calculated T1 states possessed nonemissive characters. Apart from 24 for which the nonemissive character most likely stem from the limit of the density functional employed in this work, such characteristics are in line with ref 57, in which the T1 state had been determined as a nonemissive state, and the emissive state were found, lying higher in energy. Excitation analysis and natural transition orbital (NTO) transformations were performed using TheoDORE program. 58 All calculations were carried out in C1 symmetry. 17



Page 18 of 35

Table 1: Calculated photophysical properties of the tetradentate Pt(II) complexes. kphos (104 s−1 ) Complex 14 15 16 17 18 19 20 21 22 23 25 26 27 28 29 30 31 32 33 34 35 40 41 42 43 44 45 46 47 48 49 50 51 52

Calc. 2.8 1.4 4.2 4.3 1.8 2.2 3.6 1.6 1.9 2.7 14.8 5.7 6.6 20.8 6.1 5.6 6.8 2.0 2.7 1.3 1.0 6.6 7.1 6.7 9.4 4.7 7.2 7.2 21.8 21.3 12.7 25.9 11.2 4.6

MARE 54.2c MSRE −23.5c

Expt.a 11.3 5.3 6.9 8.4 5.8 2.1 3.2 2.7 2.3 10.5 181 9.6 1.8 11.0 10.5 0.12 3.9 5.1 4.4 4.3 3.2 15.5 18.4 20.9 15.7 13.0 8.0 5.6 21.5 18.6 13.0 17.5 31.5 6.3

kisc (104 s−1 ) Calc. 4.8 18.5 6.0 3.5 1.7 3.7 1.5 1.7 13.2 2.6 8.5 8.0 8.4 3.0 8.1 38.9 12.8 1.6 13.1 3.0 5.6 2.0 2.4 3.6 3.2 3.5 3.1 3.5 7.5 8.0 12.1 5.0 6.4 2.6

Expt.a 7.5 47.4 18.7 23.9 15.8 60.4 12.7 13.1 21.5 2.3 131 3.6 11.3 6.5 7.6 118 18.9 23.4 51.1 5.7 5.5 5.7 2.0 2.3 3.9 2.1 9.0 5.8 8.8 5.2 20.3 232 18.5 1.5

ΦPL Calc. 0.37 0.07 0.41 0.55 0.52 0.38 0.70 0.48 0.12 0.51 0.63 0.42 0.44 0.87 0.43 0.13 0.34 0.56 0.17 0.31 0.15 0.77 0.75 0.65 0.75 0.57 0.70 0.67 0.74 0.73 0.51 0.84 0.64 0.64

Expt.a 0.60 0.10 0.27 0.26 0.27 0.033 0.29b , 0.20 0.17 0.097 0.82 0.58 0.74 0.14 0.63 0.58 0.001 0.17 0.26b , 0.18 0.26b , 0.08 0.43 0.37 0.73 0.90 0.90, 0.74b 0.91b , 0.80 0.86 0.47 0.49 0.85b , 0.71 0.89b , 0.78 0.83b , 0.39 0.58b , 0.07 0.97b , 0.63 0.90b , 0.81

λem (nm) Calc. 594 617 549 618 642 673 629 599 564 453 471 489 591 478 595 634 539 484 508 534 570 539 541 560 533 551 541 541 449 449 463 432 490 476

Expt.a 586 595 546 628 625 661 624 568 566 461 492 512 613 486 621 660 550 443 461 463 465 503 522 551 517 553 526 527 454 452 430 442 513 491

MAE 0.21 21.2 (15.2)d RMSE 0.23 29.3 (16.7)d MSE 0.01 8.5 (2.6)d a Experimental values are measured in an organic solvent, and are from ref 28 and therein. b Solid-state absolute PLQY, measured using an integration sphere apparatus. c Relative errors (in percentage) excluding 30 for ΦPL < 0.05, and 49 for |ΦPMMA − Φsolvent | > 0.5. d Errors excluding 32–35 based on N-heterocyclic carbenes are given in parentheses (see the text). PL PL

3.2

56.9c −29.7c

Comparison of theory and experiment

In Table 1, the calculated ΦPL and emission maximum λem of the tetradentate Pt(II) complexes and corresponding experimental data are summarized. λem has been determined from the first peak maximum of the simulated phosphorescence spectra (see Supplementary Information Figure S2; ∆E has been blue-shifted by 1100 cm−1 to account the slight underestimation of excitation energy). If available, we used the solid-state, experimental ΦPL in the comparison because undesired quenching effects such as the self-quenching and the intermolecular aggregation, which is prone for the square-planar molecular geometry, is re18


Page 19 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


duced significantly in the solid state. 10 Mean absolute errors (MAEs) are 0.21 and 21 nm for ΦPL and λem , respectively. Also, the calculated kphos and kisc have the corresponding lifetime (i.e., 1/k) in the microsecond range that is typical of phosphorescent emitters, implying the correct timescale has been obtained from the Golden Rule expressions. This is reflected in their mean absolute relative error (MARE), 54% for kphos and 57% kisc , respectively. The corresponding mean signed relative error (MSRE) are measured as −24% and −30%. In the same table, we have employed three different error metrics, MAE, root mean squared error (RMSE), and mean signed error (MSE), to capture the performance of the present method. MAE and RMSE of the calculated ΦPL are close to each other (0.21 and 0.23, respectively), which shows the absence of outliers as RMSE which, by definition, tends to penalize large errors more. By contrast, in case of λem , RMSE is 38% higher than the corresponding MAE; this is largely due to 32–35 with the common ligand scaffold consisting N-hetrocyclic carbene, showing relatively large errors (41 nm, 47 nm, 71 nm and 105 nm, respectively) owing to broad, featureless emission spectra simulated (see Figure S2), hence overestimating the λem . MAE and RMSE of λem after excluding 32–35 are brought closer to each other (15.2 nm and 16.7 nm, respectively), implying again the absence of outliers. To obtain above results, the parameter τL in the Lorentzian damping function in section 2.6 was determined by maximizing the Pearson correlation coefficient (often referred to as R-value) while the three error metrics of MAE, RMSE and MSE are kept minimized. As seen in Table 2, we obtained the best result when τL = 500 fs that corresponds to FWHM(L) of 10.6 cm−1 , all considered error metrics showed significant improvements. FWHM(L) < 100 cm−1 has been suggested and employed previously. 17,20 Table 2: Calibration of the damping parameter τL .

a

τL a (fs) Pearson MAE RMSE MSE 100 0.578 0.30 0.38 −0.29 500 0.640 0.21 0.23 0.01 1000 0.628 0.22 0.27 0.14 2500 0.574 0.28 0.36 0.27 5000 0.523 0.33 0.42 0.33 10000 0.478 0.36 0.45 0.36 Corresponding FWHM(L) are 53.1, 10.6 cm−1 , 5.3 cm−1 , 2.1 cm−1 , 1.1 cm−1 , and 0.5 cm−1 , respectively, down the column.

19



8

10

w/o mode−mixing w/ mode−mixing

0.6

107 106

0.4

105

0.2

104 3

10

0.8

Complex

RMSD (Å)

9

10 kisc (s−1)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 35

0

Figure 2: Inter-system crossing rates with and without mode-mixing effect. Complexes are arranged along the x-axis in order of increasing RMSD (gray). As mentioned in section 2.6, we have employed the displaced and distorted harmonic approximation in the calculation, and one should confirm the approximation does not alter the trend of calculations, at least in the qualitative sense. Figure 2 shows the effect of modemixing on kisc in terms of the root mean squared deviation (RMSD) of geometries, which measures the average distance between the atoms of the T1 and the S0 state geometries. The rectilinear Duschinsky matrices were utilized here, and we note that when the geometries are identical, this matrices are equivalent to the corresponding curvilinear Duschinsky matrices. 43 Because the determinants deviated from unity due to the forementioned reasons, we did not utilize curvilinear Duschinsky matrices. which were not able to be used in our calculates As clearly seen, at very small RMSD values (< 0.05 Å; the 4 left-most complexes in Figure 2) the difference between the two kisc ’s became less than 1.5%, suggesting a negligible mode-mixing effect. With moderately small RMSD values in the range of 0.05– 0.2 Å, the differences remained within an order of magnitude, apart from few cases where kisc peaked by few order of magnitudes most likely due to an ill-definition of the rectilinear Duschinsky matrix (which are found common for RMSD greater than 0.2 Å). Since no obvious outliers of kisc have been observed with this approximation, the mode-mixing effect on kisc for the tetradentate Pt(II) complexes can be assumed relatively small, suggesting that the de-excitation dynamics of the considered Pt(II) complexes can be effectively described within the displaced and distorted harmonic oscillator model.

20


Page 21 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3: Molecules whose nonradiative characteristics are investigated.

3.3

Nonradiative characteristics of Pt(Ph2 N2 O2 ) and Pt(tBu2 N2 O2 )

To shed light on the source of rapid nonradiative T1

S0 ISC in some phosphorescent emit-

ters, the complexes 14 [Pt(Ph2 N2 O2 )] and 15 [Pt(tBu2 N2 O2 )] are selected for detailed analyses (Figure 3). 59 As evident from Table 1, kisc has increased by approximately 4-fold for Pt(tBu2 N2 O2 ) (18.5 × 104 s−1 ) from Pt(Ph2 N2 O2 ) (4.8 × 104 s−1 ) when the bridge between the pyridinyl rings is removed. In other words, the ISC of complexes based on the common OˆNˆNˆO ligand scaffold can be suppressed significantly when the two pyridinyl rings are enforced by a rigid fused aromatic bridge. The effect of above structural modification is well portrayed in Figure 4 where the nonBorn–Oppenheimer coupling factors Tκκ and the reorganization energies λκ ≡ 21 d2κ ωκ2 as well as the diagonal ISC rates [kisc ]κκ have been decomposed for each normal mode. It is clearly seen that Tκκ and λκ for Pt(Ph2 Bu2 N2 O2 ) are much smaller than that of Pt(tBu2 N2 O2 ) by at least an order, and the normal mode ν157 of the latter at 1587 cm−1 , shown in Figure 5 similar to the C–C inter-ring vibrations between the two pyridinyl rings, is characterized by having the maximum values in both fields: T157,157 = 679 cm−2 and λ157 = 732 cm−1 . This indicates that the larger degree of geometrical displacements along ν157 between the states coincides with a higher non-Born-Oppenheimer coupling factor. In the Golden Rule expression, this leads to a fast nonradiative transition between the potential surfaces since the curve-crossing point is effectively lowered. There are, of course, other non-obvious, latent factors affecting the final rate kisc since the contribution of ν157 is merely 39.2% (Figure 4b). Even though being a major component, this does not explain the overall behavior by itself and one should explicitly solve eq 35 for kisc . We note that the diagonal contributions to kisc 21



0

Pt(tBu2N2O2)

800

400 800

ν157

400 0

0

500

b)

1000 1500 0 500 1000 Normal mode frequency (cm−1)

λκ (cm−1)

Pt(Ph2N2O2)

[kisc]κκ (s−1)

a) Tκκ (cm−2)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

80000 60000

Pt(Ph2N2O2)

Pt(tBu2N2O2)

diagonal = 98.4%

diagonal = 96.7%

ν157

40000 20000 0

1500

Page 22 of 35

0

500

500 1000 1000 1500 0 Normal mode frequency (cm−1)

1500

Figure 4: (a) Electronic non-Born–Oppenheimer electronic coupling factors Tκκ and reorganization energies λκ , and (b) diagonal contributions [kisc ]κκ to the final kisc for each normal mode. Negligible values are obtained for the terminal hydrogen stretches at ∼3000 cm−1 and above, thus omitted for brevity.

Figure 5: Illustration of the normal mode ν157 at 1587 cm−1 in Pt(tBu2 N2 O2 ) with atomic numbering scheme. Hydrogens are omitted. (i.e.,

P

κ [kisc ]κκ /kisc )

for both cases are over 95% as expected for the displaced and distorted

harmonic approximation used. Table 3 summarizes the large changes in bond length between the S0 to the T1 states in Pt(tBu2 N2 O2 ). Most bond contractions and elongations are concentrated in the pyridinyl rings, and, although not shown in the table, dihedral torsions less than 10 degrees are additionally detected along the Pt–O bonds. The largest bond contraction (approximately 0.06 Å) occurs between C5 –C8 upon excitation to the T1 state; this explains the largest reorganization energy along the mode ν157 . Furthermore, the analysis on the electronic structure of the excited states employing NTO representation showed that the T1 state can be characterized by a mixed metal-to-ligand (MLCT) and inter-ligand charge-transfer (ILCT) transitions, populating the π ∗ -orbital across the pyridinyl rings (Figure 6) of local πbonding characters along the C5 –C8 bond. The NTO pair of an electronic state illustrates the

22


Page 23 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 3: Changes in bond length (Å) between the S0 and the T1 state in Pt(tBu2 N2 O2 )

a

Bonda S0 state T1 state Differenceb N4 –C5 1.360 1.391 0.030 C5 –C6 1.386 1.416 0.030 C8 –C9 1.386 1.416 0.030 C8 –N10 1.360 1.391 0.030 C2 –C7 1.392 1.417 0.024 C12 –C13 1.392 1.417 0.024 C6 –C7 1.398 1.378 −0.020 C9 –C13 1.398 1.378 −0.020 C5 –C8 1.471 1.418 −0.053 Indices follow from Figure 5. b Differences greater than 0.02 Å are listed.

T1 (S1 ) : 0.99(0.99)

S2 : 0.97

LUNTO

S3 : 0.97

S7 : 0.84

S12 : 0.85

Figure 6: Natural transition orbitals of the excited states of Pt(tBu2 N2 O2 ). HONTOs depicted in red/blue share the same LUNTO composition (yellow/green). Corresponding eigenvalues are given. excitation in a simple orbital-like picture: highest-occupied (HO)NTO and lowest-unoccupied (LU)NTO describe the hole and the excited electron state, respectively. 60 Next, to understand the nature of the higher non-Born–Oppenheimer coupling factor of ν157 , we examined the convergence of the vector t157 with respect to the included excited singlet states from the S1 through the S100 state in the perturbative expansion eq 30 (Figure 7). With first 50 singlet excited states, the convergence within 1 cm−1 has reached, attesting our choice of the singlet manifold truncation in this work. The first-order electronic nonadiabatic coupling in eq 30 can be recasted using perturbation theory: 61 hSn | ∂Q∂ κ Vˆne |S0 i ∂ hSn | S0 i = , ∂Qκ ES0 − ESn

(37)

where Vˆne denotes the nucleus-electron attraction operator. The condition ET1 > 2 eV for

23



S12

30 |t157| (cm−1)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

S50 S100

S7 S2

20 10 0 0

1

2 3 ES − ET (eV) n

4

5

1

Figure 7: Convergence profile of |t157 | with respect to the number of excited singlet states included in a sum-over-states calculation. phosphorescent emitters in the visible photons, the denominator in eq 30 can be estimated as (∆2 + 2∆) eV where ∆ = ESn − ET1 . Assuming that the numerator is kept similar order, for instance, the contribution from Sn with ∆ = 0.5 eV is roughly 20-fold larger than one with ∆ = 4 eV. Therefore, the expansion would converge roughly quadratically, and the elements in the coupling matrix T quartically. The contributions from the nearby excited singlet states are, therefore, of essential importance in the nonradiative dynamics to the S0 state. Table 4 compiles the calculated spin–orbit and nonadiabatic couplings in Pt(tBu2 N2 O2 ). The S2 state in Figure 7 has been observed as the major contributor to t157 , along with minor contributions from the S7 and the S12 states, and its largest atomic norms reside on the bridging atoms, C5 and C8 . Moreover, the corresponding spin–orbit coupling with the T1 state is also relatively strong (reaching 338 cm−1 ) due to fact that the dxz and dyz orbitals are involved in the T1 and the S2 state, respectively, (Figure 6) with a relatively small energy difference of 0.76 eV. These factors allow an effective mixing of the S2 state into the spin-mixed solution of the T1 state, leading to the rapid nonradiative relaxation along the inter-ring vibrational coordinates. The S7 and the S12 states demonstrate similar characteristics, however, the larger energy gaps of 1.6 eV and 2.2 eV between these states, respectively, have reduced their contributions. The S1 state, on the contrary, is associated with the same Pt d-orbital as the T1 state, leading to negligible spin–orbit coupling and subsequently, a little participation in the nonradiative dynamics. Lastly, HONTO of the S3

24


Page 24 of 35

Page 25 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 4: Spin–orbit coupling (cm−1 ) and norm of the nonadiabatic coupling vector (bohr−1 ) for Pt(tBu2 N2 O2 ). The electronic energy difference ESn − ET1 is given in eV units. n | hS0 | ∂S i |a ∂ξ 0.519 0.576 0.250 0.564

n 1 2 3 4

ESn − ET1 0.27 0.76 1.35 1.38

ˆ SO |T1 i | | hSn |H 0.9 337.6 2246.4 368.9

5

1.38

125.9

0.635

6

1.49

14.4

0.555

7

1.64

503.9

0.646

8 9 10 11 12 13 a

1.77 182.2 0.269 1.85 26.0 0.420 2.00 109.4 0.426 2.00 61.9 0.426 2.23 422.0 0.495 2.47 253.3 0.521 Calculated over 3N Cartesian components. b

Norms of each elementb N4 , N10 (0.20); C3 , C11 , C6 , C9 (0.15); C7 , C13 (0.14); Pt1 (0.11) C5 , C8 (0.26); C13 , C7 , N10 (0.14); N4 (0.13); C9 , C6 (0.12); C12 , C2 (0.11) Pt1 (0.17); N10 , N4 (0.10) C18 (0.21); C16 (0.19); O27 (0.18); C20 , C21 (0.17); C17 , C19 (0.16); C3 (0.15); Pt1 (0.13); C2 (0.10) C25 (0.23); C15 (0.21); C11 , C23 (0.18); C24 , C26 , C22 (0.17); O28 (0.16); C12 (0.14); C18 (0.12); C16 , C3 (0.11) C3 , C11 (0.17); C6 , C2 , C12 :, C9 (0.16); C5 , C8 , Pt1 , N4 (0.13); N10 (0.12); C7 , C13 (0.11) N4 , N10 (0.20); C3 , C11 (0.19); C7 , C13 (0.18); C6 , C9 , C2 (0.17); C12 (0.16); C5 (0.15); C8 (0.14) C3 , C11 (0.11) C7 , C13 (0.14); C3 , C11 , C6 (0.13); C9 , N4 , N10 (0.12); C2 , C12 (0.10) N4 (0.18); C18 (0.16); C17 (0.15); C16 , C3 , C13 (0.14); C20 (0.13); C19 (0.12) N10 (0.18); C25 (0.16); C15 (0.15); C26 , C11 (0.14); C23 (0.13); C24 (0.11) C5 , C8 (0.22); C9 , C6 (0.16); C13 , C7 (0.15) N4 , N10 , C7 , C8 , C5 , C6 , C2 (0.13); C9 , C13 , C12 (0.12) Indices follow from Figure 5. Norms greater than 0.1 bohr−1 are given in parentheses.

state (Figure 6) shows an association with the dz2 orbital, leading to the largest spin–orbit coupling of 2246 cm−1 . As the nonadiabatic coupling is almost localized on the motion of Pt and the coordinating N-atoms, the S3 state also has a negligible contribution to the ISC transition along the ν157 coordinates. In summary, the S2 state in Pt(tBu2 N2 O2 ) plays a pivotal role for the relatively large kisc , and the complicated interplay among small ∆E, large spin–orbit coupling with the T1 , and appreciable nonadiabatic coupling with the S0 state, has all led to the faster ISC along the C–C inter-ring vibration, which also carry large portion of the reorganization energy between the T1 state and the S0 state. An introduction of a rigid aromatic bridge in Pt(Ph2 N2 O2 ) has extended the conjugation plane on the ligand scaffold, and has led to more delocalized character of LUNTO and to increased structural stiffness upon the excitation, resulting much smaller nonadiabatic coupling, and consequently slow nonradiative relaxation.

4

Conclusion

In this work, we have devised a quantitative prediction method for the phosphorescence quantum yield by employing the correlation function approach. Specifically, the two main 25



components of the PLQY, namely, the rates of phosphorescence kphos and the inter-system crossing kisc from the excited T1 state to the ground S0 state, are obtained through the correlation functions calculated in time domain by Fourier transformation of the Golden Rule rate equations. The spin–orbit coupling is efficiently treated in the first-order perturbation theory using the relativistic effective core potential in order to generate spin–orbit corrected, two-component T1 degenerate substates within the single-excitation theory. The electronic transition dipole moments and the non-Born-Oppenheimer coupling matrix elements between the S0 state and the T1 substates are computed, which are utilized in calculations of kphos and kisc , respectively. Using the substate-averaging at high-temperature limit and the timereversal symmetry constraints, we show that the ISC and the IC processes can be considered identically in the perspective of implementation by eliminating complex-valued modifications in the algorithm. Results have been obtained for the photophysical properties of 34 Pt(II) complexes designed for the OLED applications, and have shown good agreement with experiments by producing consistent predictions over diverse ligand scaffolds with an average error of 0.21. Further analysis based on the decomposition of contributions from the nuclear vibrations and from the excited-state electronic structure to the non-Born–Oppenheimer coupling has also been presented as a way to deeper insight into the source of rapid nonradiative characteristics. Securing structural rigidity in the chelating tetradentate ligand by an aromatic bridge are shown as an example of chemical modification towards higher PLQY of the Pt(II)-based OLED emitters via reducing the nonradiative decay channel facilitated by an inter-ring vibration. As a final remark, we underline that our model can be extended to the case of ISC from the triplet metal-centered state, 25 although, in this case, the accuracy of our model relies severely on the applicability of harmonic and Condon approximations since the metalcentered state usually results in large structural changes. Our implementation provides additional means of narrowing down optimal emitters in an a priori fashion in the development of novel phosphorescent emitters for OLED device,

26


Page 26 of 35

Page 27 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


as well as presenting a new platform to design highly efficient phosphorescent emitters with desired photophysical properties for industrial applications.

Acknowledgement The authors thank Professors Yoon Sup Lee and Young-Min Rhee (KAIST) and Dongwook Kim (Kyonggi University) for helpful discussions. Computational resources were provided by the supercomputing center of SAIT.

Supporting Information Available Chemical structures, optimized geometries, vibrational frequencies, phosphorescence spectra, and non-Born–Oppenheimer coupling factors of all the 34 structures used in this work. These materials are available free of charge via the Internet at http://pubs.acs.org.

References (1) Yersin, H. Highly Efficient OLEDs with Phosphorescent Materials; Wiley, 2008. (2) Baldo, M. A.; O’Brien, D. F.; You, Y.; Shoustikov, A.; Sibley, S.; Thompson, M. E.; Forrest, S. R. Highly Efficient Phosphorescent Emission From Organic Electroluminescent Devices. Nature 1998, 395, 151–154. (3) Tsujimura, T. OLED Displays: Fundamentals and Applications; Oxford Woodhead Publishing, 2012. (4) Buckley, A. Organic Light-Emitting Diodes (OLEDs): Materials, Devices and Applications; Wiley, 2013. (5) Xu, H.; Chen, R.; Sun, Q.; Lai, W.; Su, Q.; Huang, W.; Liu, X. Recent Progress in

27



Metal–Organic Complexes for Optoelectronic Applications. Chem. Soc. Rev. 2014, 43, 3259–3302. (6) Chi, Y.; Chou, P.-T. Transition-Metal Phosphors With Cyclometalating Ligands: Fundamentals and Applications. Chem. Soc. Rev. 2010, 39, 638–655. (7) Brütting, W.; Adachi, C. Physics of Organic Semiconductors; Wiley-VCH, 2013. (8) Minaev, B.; Baryshnikov, G.; Agren, H. Principles of Phosphorescent Organic Light Emitting Devices. Phys. Chem. Chem. Phys. 2014, 16, 1719–1758. (9) Penfold, T. J.; Gindensperger, E.; Daniel, C.; Marian, C. M. Spin-Vibronic Mechanism for Intersystem Crossing. Chem. Rev. 2018, 118, 6975–7025. (10) Ma, B.; Djurovich, P. I.; Thompson, M. E. Excimer and Electron Transfer Quenching Studies of a Cyclometalated Platinum Complex. Coord. Chem. Rev. 2005, 249, 1501– 1510. (11) Baryshnikov, G.; Minaev, B.; Ågren, H. Theory and Calculation of the Phosphorescence Phenomenon. Chem. Rev. 2017, 117, 6500–6537. (12) Ianconescu, R.; Pollak, E. Photoinduced Cooling of Polyatomic Molecules in an Electronically Excited State in the Presence of Dushinskii Rotations. J. Phys. Chem. A 2004, 108, 7778–7784. (13) Niu, Y.; Peng, Q.; Deng, C.; Gao, X.; Shuai, Z. Theory of Excited State Decays and Optical Spectra: Application to Polyatomic Molecules. J. Phys. Chem. A 2010, 114, 7817–7831. (14) Ito, A.; Meyer, T. J. The Golden Rule. Application for Fun and Profit in Electron Transfer, Energy Transfer, and Excited-State Decay. Phys. Chem. Chem. Phys. 2012, 14, 13731.

28


Page 28 of 35

Page 29 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(15) Etinski, M.; Tatchen, J.; Marian, C. M. Time-Dependent Approaches for the Calculation of Intersystem Crossing Rates. J. Chem. Phys. 2011, 134, 154105. (16) Kleinschmidt, M.; van Wüllen, C.; Marian, C. M. Intersystem-Crossing and Phosphorescence Rates in fac-IrIII (ppy)3 : A Theoretical Study Involving Multi-Reference Configuration Interaction Wavefunctions. J. Chem. Phys. 2015, 142, 094301. (17) Peng, Q.; Niu, Y.; Shi, Q.; Gao, X.; Shuai, Z. Correlation Function Formalism for Triplet Excited State Decay: Combined Spin–Orbit and Nonadiabatic Couplings. J. Chem. Theory Comput. 2013, 9, 1132–1143. (18) Baiardi, A.; Bloino, J.; Barone, V. General Time Dependent Approach to Vibronic Spectroscopy Including Franck–Condon, Herzberg–Teller, and Duschinsky Effects. J. Chem. Theory Comput. 2013, 9, 4097–4115. (19) Tapavicza, E.; Furche, F.; Sundholm, D. Importance of Vibronic Effects in the UV–Vis Spectrum of the 7,7,8,8-Tetracyanoquinodimethane Anion. J. Chem. Theory Comput. 2016, 12, 5058–5066. (20) de Souza, B.; Neese, F.; Izsák, R. On the Theoretical Prediction of Fluorescence Rates From First Principles Using the Path Integral Approach. J. Chem. Phys. 2018, 148, 034104. (21) de Souza, B.; Farias, G.; Neese, F.; Izsák, R. Predicting Phosphorescence Rates of Light Organic Molecules Using Time-Dependent Density Functional Theory and the Path Integral Approach to Dynamics. J. Chem. Theory Comput. 2019, 15, 1896–1904. (22) Gourlaouen, C.; Eng, J.; Otsuka, M.; Gindensperger, E.; Daniel, C. Quantum Chemical Interpretation of Ultrafast Luminescence Decay and Intersystem Crossings in Rhenium(I) Carbonyl Bipyridine Complexes. J. Chem. Theory Comput. 2015, 11, 99–110.

29



(23) Eng, J.; Gourlaouen, C.; Gindensperger, E.; Daniel, C. Spin-Vibronic Quantum Dynamics for Ultrafast Excited-State Processes. Acc. Chem. Res. 2015, 48, 809–817. (24) Escudero, D. Quantitative Prediction of Photoluminescence Quantum Yields of Phosphors From First Principles. Chem. Sci. 2016, 7, 1262–1267. (25) Sajoto, T.; Djurovich, P. I.; Tamayo, A. B.; Oxgaard, J.; Goddard, W. A.; Thompson, M. E. Temperature Dependence of Blue Phosphorescent Cyclometalated Ir(III) Complexes. J. Am. Chem. Soc. 2009, 131, 9813–9822. (26) Zhang, X.; Jacquemin, D.; Peng, Q.; Shuai, Z.; Escudero, D. General Approach to Compute Phosphorescent OLED Efficiency. J. Phys. Chem. C 2018, 122, 6340–6347. (27) Dolg, M.; Cao, X. Relativistic Pseudopotentials: Their Development and Scope of Applications. Chem. Rev. 2012, 112, 403–480. (28) Li, K.; Ming Tong, G. S.; Wan, Q.; Cheng, G.; Tong, W.-Y.; Ang, W.-H.; Kwong, W.L.; Che, C.-M. Highly Phosphorescent Platinum(II) Emitters: Photophysics, Materials and Biological Applications. Chem. Sci. 2016, 7, 1653–1673. (29) Hofbeck, T.; Yersin, H. The Triplet State of fac-Ir(ppy)3 . Inorg. Chem. 2010, 49, 9290– 9299. (30) Runge, E.; Gross, E. K. U. Density-Functional Theory for Time-Dependent Systems. Phys. Rev. Lett. 1984, 52, 997–1000. (31) Gross, E. K. U.; Kohn, W. Local Density-Functional Theory of Frequency-Dependent Linear Response. Phys. Rev. Lett. 1985, 55, 2850–2852. (32) Hirata, S.; Head-Gordon, M. Time-Dependent Density Functional Theory Within the Tamm–Dancoff Approximation. Chem. Phys. Lett. 1999, 314, 291–299. (33) Daniel, C. Photochemistry and Photophysics of Transition Metal Complexes: Quantum Chemistry. Coord. Chem. Rev. 2015, 282-283, 19–32. 30


Page 30 of 35

Page 31 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(34) McClure, D. S. Spin-Orbit Interaction in Aromatic Molecules. J. Chem. Phys. 1952, 20, 682–686. (35) Powell, B. Theories of Phosphorescence in Organo-Transition Metal Complexes – From Relativistic Effects to Simple Models and Design Principles for Organic Light-Emitting Diodes. Coord. Chem. Rev. 2015, 295, 46–79. (36) Dobbs, K. D.; Feldman, J.; Marshall, W. J.; McLain, S. J.; McLaren, C. D.; Meth, J. S.; Vo, G. D.; Wang, Y. Phosphorescent Iridium(III) Complexes of Cyclometalated 5Aryl-1H-1,2,4-Triazole Ligands: Structural, Computational, Spectroscopic, and Device Studies. J. Phys. Chem. C 2014, 118, 27763–27771. (37) Tunell, I.; Rinkevicius, Z.; Vahtras, O.; Sałek, P.; Helgaker, T.; Ågren, H. Density Functional Theory of Nonlinear Triplet Response Properties With Applications to Phosphorescence. J. Chem. Phys. 2003, 119, 11024–11034. (38) Parker, S. M.; Rappoport, D.; Furche, F. Quadratic Response Properties from TDDFT: Trials and Tribulations. J. Chem. Theory Comput. 2018, 14, 807–819. (39) Kramers, H. A. Théorie Générale De La Rotation Paramagnétique Dans Les Cristaux. Proc. R. Acad. Sci. Amsterdam 1930, 33, 959–972. (40) Rösch, N. Time-Reversal Symmetry, Kramers’ Degeneracy and the Algebraic Eigenvalue Problem. Chem. Phys. 1983, 80, 1–5. (41) May, V.; Kühn, O. Charge and Energy Transfer Dynamics in Molecular Systems; WileyVCH, 2011. (42) Piela, L. Ideas of Quantum Chemistry; Elsevier B.V., 2014. (43) Reimers, J. R. A Practical Method for the Use of Curvilinear Coordinates in Calculations of Normal-Mode-Projected Displacements and Duschinsky Rotation Matrices for Large Molecules. J. Chem. Phys. 2001, 115, 9103–9109. 31



(44) Chang, J.-L. A New Formula to Calculate Franck–Condon Factors for Displaced and Distorted Harmonic Oscillators. J. Mol. Spectrosc. 2005, 232, 102–104. (45) Frigo, M.; Johnson, S. G. The Design and Implementation of FFTW3. Proc. IEEE 2005, 93, 216–231 (46) Kalinowski, J.; Fattori, V.; Cocchi, M.; Williams, J. G. Light-Emitting Devices Based on Organometallic Platinum Complexes as Emitters. Coord. Chem. Rev. 2011, 255, 2401–2425. (47) Börgel, J.; Campbell, M. G.; Ritter, T. Transition Metal d-Orbital Splitting Diagrams: An Updated Educational Resource for Square Planar Transition Metal Complexes. J. Chem. Educ. 2016, 93, 118–121. (48) Gourlaouen, C.; Daniel, C. Spin–Orbit Effects in Square-Planar Pt(II) Complexes With Bidentate and Terdentate Ligands: Theoretical Absorption/Emission Spectroscopy. Dalt. Trans. 2014, 43, 17806–17819. (49) Furche, F.; Ahlrichs, R.; Hättig, C.; Klopper, W.; Sierka, M.; Weigend, F. Turbomole. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2014, 4, 91–100. (50) Perdew, J. P.; Ernzerhof, M.; Burke, K. Rationale for Mixing Exact Exchange With Density Functional Approximations. J. Chem. Phys. 1996, 105, 9982–9985. (51) Adamo, C.; Barone, V. Toward Reliable Density Functional Methods Without Adjustable Parameters: The PBE0 Model. J. Chem. Phys. 1999, 110, 6158–6170. (52) Sierka, M.; Hogekamp, A.; Ahlrichs, R. Fast Evaluation of the Coulomb Potential for Electron Densities Using Multipole Accelerated Resolution of Identity Approximation. J. Chem. Phys. 2003, 118, 9136–9148. (53) Weigend, F.; Ahlrichs, R. Balanced Basis Sets of Split Valence, Triple Zeta Valence

32


Page 32 of 35

Page 33 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297. (54) Figgen, D.; Peterson, K. A.; Dolg, M.; Stoll, H. Energy-Consistent Pseudopotentials and Correlation Consistent Basis Sets for the 5d Elements Hf–Pt. J. Chem. Phys. 2009, 130, 164108. (55) Send, R.; Furche, F. First-Order Nonadiabatic Couplings From Time-Dependent Hybrid Density Functional Response Theory: Consistent Formalism, Implementation, and Performance. J. Chem. Phys. 2010, 132, 044107. (56) Klamt, A. The COSMO and COSMO-RS Solvation Models. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2018, 8, e1338. (57) Kui, S. C.; Chow, P. K.; Tong, G. S. M.; Lai, S.-L.; Cheng, G.; Kwok, C.-C.; Low, K.H.; Ko, M. Y.; Che, C.-M. Robust Phosphorescent Platinum(II) Complexes Containing Tetradentate OˆNˆCˆN Ligands: Excimeric Excited State and Application in Organic White-Light-Emitting Diodes. Chem. Eur. J. 2013, 19, 69–73. (58) Plasser, F.; Wormit, M.; Dreuw, A.; Bäppler, S. A.; Wormit, M.; Dreuw, A. New Tools for the Systematic Analysis and Visualization of Electronic Excitations. I. Formalism. J. Chem. Phys. 2014, 141, 024107. (59) Lin, Y.-Y.; Chan, S.-C.; Chan, M. C. W.; Hou, Y.-J.; Zhu, N.; Che, C.-M.; Liu, Y.; Wang, Y. Structural, Photophysical, and Electrophosphorescent Properties of Platinum(II) Complexes Supported by Tetradentate N2 O2 Chelates. Chem. Eur. J. 2003, 9, 1263–1272. (60) Martin, R. L. Natural Transition Orbitals. J. Chem. Phys. 2003, 118, 4775. (61) Chernyak, V.; Mukamel, S. Density-Matrix Representation of Nonadiabatic Couplings

33



in Time-Dependent Density Functional (TDDFT) Theories. J. Chem. Phys. 2000, 112, 3572–3579.

34


Page 34 of 35

Page 35 of 35

Graphical TOC Entry

kisc

kphos 0

10000

Phosphorescence Quantum Yield

Inter-system crossing

0 Energy (cm−1)

10000

35

Log(ISC rate)

Phosphorescence Photon emission rate

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ΦPL =

20000


kphos kphos + kisc

Predicting Phosphorescence Quantum Yield for Pt(II)-Based OLED

Recommend Documents