The newest display technologies could not be realized without the developments in organic light-emitting diodes (OLEDs). Despite their commercial success, there are still issues regarding the efficiency and lifetime of these devices to bring the technology to the next level. Similar (or more drastic issues) need to be resolved to develop a LED technology based on perovskite light-emitting materials. It is especially the blue LED that is lagging behind the performance and reliability of the green and red pixels. This affects the overall display operation, pixel layout choice, and performance.¹ Therefore, it is important to better understand the device physics, in general, and the origins of degradation in OLEDs. Novel insight into the operating mechanisms will allow for the design and selection of new materials and provide guidance to optimize the device structure.

Besides the traditional solution or film characterization techniques, such as cyclic voltammetry^2–4and ultraviolet photoelectron spectroscopy (UPS),⁵ new materials need to be thoroughly characterized in full multilayer LEDs or specific layer stacks to assess their suitability in real applications. This is also important as material properties may depend on the device fabrication details as well as the specific layer sequence.^6–10 Thus, reliable material characterization in organic and perovskite LEDs is a key requirement. Yet, when devices are used as characterization test-beds, the material parameters are still extracted from single measurements on specific samples. This approach is time-consuming, requires high material consumption and the determined parameters are of limited validity because they might depend on the explicit device structure or on the peculiar analytical model that was selected to fit the data.^11–13

The preferred way to determine material parameters in organic and perovskite LEDs should be, therefore, to combine numerical simulations with experiments. Often, these comparisons were solely focusing on steady-state analysis.^14–16 In such cases, it is necessary to be careful that the analyzed material parameters are unique.^17–19 In order to increase the reliability of device simulation, it is inevitable to include complementary measurement techniques in time and/or frequency domain.^19–21

where d is the thickness of the device and Vbi the built-in potential. Usually the transit time is approximated by the delay time td between the voltage turn-on and the light onset which might lead to an overestimation of the mobility. ³⁶ This simple formula was developed for single-layer devices which are not injection limited and neglects the spatial variation of the electric field due to the accumulated charge. ^36–38 In multi-layer bipolar devices, the relevant thickness parameter d corresponds to the distance which the slower charge carrier needs to travel in order to meet and recombine with the counter-charge. Even though in state-of-the-art OLEDs the onset analysis only yields an apparent mobility, the technique can nevertheless be very useful to qualitatively understand basic charge transport phenomena. ³⁹ The analysis of the TEL rise dynamics is not discussed here in detail for sake of brevity, even though fast and subsequent slow rise can be assigned to the two carriers in the emission layer36 and sometimes a TEL onset overshoot can be observed due to rapid depletion of internally accumulated charges. ⁴⁰

The TEL decay, however, is generally related to exciton dynamics and interpreted similarly to transient photoluminescence (TRPL). A single exponential decay fit yields an exciton lifetime that can be related to the radiative and non-radiative recombination rates. In contrast to TRPL also charge carriers are present in a TEL experiment. Therefore the decay can be affected by exciton-polaron interactions41 and delayed recombination of charge carriers. ²⁰ Such TEL data have been analyzed in depth do understand different losses in OLED devices ^20,42,43 as well as changes during device degradation. ⁴⁴

Parameter that can be extracted from a TE analyis:
Luminescence Lifetime; Field Dependent Charge Carrier Mobility.

In all CELIV techniques, the characteristic current peak area, namely the time integral of the current overshoot, can be analyzed to estimate the amount of accumulated charge in the device. The peak time t _max is related to the charge carrier mobility µ following the analytically derived formula

where d is the thickness of the device and A is the employed voltage ramp. Several modifications of this formula have been presented to account for commonly neglected effects such as the non-uniform field and charge distribution, series resistance and recombination. ^51–56 As for the TEL rise, the use of multilayer OLED structures complicates or even excludes the application of the analytical formula to extract charge carrier mobilities accurately and assign the mobility value to a specific material or carrier type. In such devices the technique should rather be used as qualitative assessment of charge carrier dynamics, e.g.to analyze differences between devices or to identify systematic changes upon degradation.

Because for state-of-the-art multilayer OLED stacks the above presented analytical formulas are only of limited validity to determine material parameters, full device modelling should be used preferably. The electrical simulations applied in the present study are based on the drift-diffusion approach which has been presented in many previous publications. ^17,31,38

The equations can be solved in steady-state, transient and frequency domain to reproduce the measured curves of the electrical experiments. For electro-optical simulations the calculated exciton densities are transferred to the optical solver which is based on the dipole emission and transfer matrix models. ^57,58 The commercially available simulation software Setfos is used in this study. ⁵⁹

The required inputs for the simulation are the layer sequence with respective thicknesses, the refractive indices, and the electrical and excitonic material parameters of each layer. Here, for simplification, and motivated by the band diagram in Figure 1(a), we choose discrete energy levels for transport and traps. For the emitter layer, charge transport is assumed to occur on the mCBP HOMO (holes) and the LUMO of the emitter (electrons). In this example the OLED stack consist mostly of standard materials without any electrical doping.

We set up a protocol, which is shown as flow chart in Figure 2(a). This procedure shall be employed to simulate the measured data on the OLED devices presented in Figure 1.

A full electrical characterization of an OLED typically consists of various steady-state, transient and impedance measurements. As a first step, a selection of relevant results has to be performed. In general, it is recommended to consider as many experimental techniques as possible. One should thus try to select measurements that yield complementary data, but avoid to include multiple similar experiments that contain the same information. For instance, injection-CELIV data at one or two different offset voltages is usually enough. In a second step, we determine those device and material parameters that can be extracted from simple post-processing techniques (e.g. series resistance). Next, the initial guess simulation is set up using measured, assumed, or literature values for the missing material parameters. Importantly, in order to minimize the number of free parameters, and in agreement with Occam’s razor, the initial model should be as simple as possible. Applied to the present case study, this means for example that we shall neglect any trap states at first, and only introduce them in case the simple model fails to properly describe the measured data. After a first simulation of the initial situation one should vary all fitting parameters to learn about their effect on the different measurement techniques. This also allows to manually optimize these parameters in a next step in order to match the experimental data as well as possible. If the agreement after the optimization is not sufficient, new model parameters should be introduced, for example trap states, and the optimization should be done again. Finally, once the agreement between simulation and measurements is satisfying, the fit quality should be analyzed. This task consists of a parameter correlation ^17,18 and sensitivity ^19,60 analysis.

We will guide the reader through the workflow illustrated in Figure 2 on the basis of real measurements on OLEDs. The selected experimental data which we want to reproduce with device simulations is shown in Figure 3.

The IV curve shows an onset at about 2.2 V for all devices with 40 nm thick NPB. As this parameter is mostly determined by the built-in potential it could be speculated that the effect is caused by the fabrication process. Furthermore, the current in the on state (above 3 V) is systematically reduced with increasing NPB thickness. This indicates that the current at high voltages is limited by hole transport inside the NPB layer.

The capacitance and conductance at V _DC = 3 V are shown in  Figure 3 (b) and (c), respectively. In the frequency range between 10 ⁵ – 10 ⁶ Hz the capacitance reflects the geometrical capacitance of the full device which is also determined at V _DC = 0 V. Towards lower frequencies, a capacitance rise is observed which is related to injection of charge carriers. This characteristic frequency f _on (see Figure 3b) is independent of the NBPhen layer thickness but shows a dependence on NPB thickness, indicating that the onset is related to hole injection and transport. A second rise is observed at a frequency around 1 kHz for all devices and needs to be linked to a process with slower dynamics.

The conductance between 10 ⁴ – 10 ⁵ Hz shows a similar thickness dependence as the IV characteristics at higher voltages. Towards lower frequencies it converges to one conductance value for all devices that have the same NBPhen thickness which is in agreement with the current in the IV measurement at the same voltage (3 V). We attribute this behavior to charge carriers (holes) that are already present in the device at 3 V and accumulate either at the TCTA/EML or at the EML/NBPhen interface. The remaining NBPhen (and EML) layer(s) is (are) thus limiting the conductivity at this voltage.

Injection-CELIV was performed at several offset voltages. For discussion and simulation purposes we chose an offset voltage of 3 V which is just above the turn-on voltage of the OLED. At this low offset voltage only small amounts of charges are injected that accumulate at some internal interface and are then extracted by the reverse voltage ramp. As outlined in section II, the peak time is related to the charge carrier mobility. For devices with the same NPB but different NBPhen thickness, the peak occurs at  the same time (Figure 3(b)). Only the peak area changes which can be explained by the different amount of charge that was injected. The NPB thickness variation results in a clear peak shift. Judging from the observed independence of t _max from the electron transport layer thickness, we can already attribute the peak to an effective hole transport time. However, without knowing the interface at which charges accumulate prior to the CELIV ramp we cannot estimate a reliable value for the hole mobility. Additionally, for multilayer structures one can only obtain an average hole mobility which is not directly related to one specific layer. It is worth noting that the measured displacement current plateau value for times t > 15 μs is in perfect agreement with the geometric capacitance of the devices measured by impedance spectroscopy.

The TEL rise and decay for a voltage pulse of 7 V are shown in  Figure 3 (e) and (f), respectively. The TEL turn-on time is related to charge transport and injection and shows a similar NPB thickness dependence as the IV curve, the first C-f onset and the injection-CELIV peak time. As for the injection-CELIV peak we cannot directly extract a mobility using the analytical formula. However, from the thickness trend we can see that the TEL onset is limited by hole transport to the emission zone. This qualitative understanding can already provide first insight in view of OLED device design and optimization. In addition to the onset time trends, a small luminance overshoot is observed for all devices before they reach their steady-state luminance at the end of the pulse. Below we will illustrate that this behavior can be related to a trapping effect in one of the hole transport layers.

In the TEL decay we do not see any difference among the devices of varying thickness. By fitting a single exponential decay to the first decay we can extract a radiative lifetime of about 5 us. This agrees well with the delayed fluorescence lifetime of the investigated TADF emitter measured by transient photoluminescence (PL) of the pure emitter film (data not shown). Interestingly, the modified cavity defined by the different layer thicknesses does not influence this decay lifetime. This finding might be related to the TADF emitter for which the reverse intersystem crossing is dominating the decay dynamics.

The above qualitative discussion demonstrates how much more information can be obtained for a specific OLED stack if consistent data with varied layer thickness is available. When fitting the experimental data in the following, the value of the complementary data set becomes even more evident.

As indicated in Figure 2(a) some device and material parameters can directly be extracted from the measurements and should be determined before setting up the simulation. The first parameters are the external series resistance Rs and parallel resistance Rp of the cells. As they depend on the contact quality as well as some shunts that may occur for one or the other cell, these parameters are extracted for each cell and not considered as global parameters. Nevertheless, it turns out that these values are very similar for all the tested devices reported here. R _s is determined by fitting the C-f measurement at V _DC = 0 V with an equivalent circuit model (Rs+[Rp||C]). A single diode model was adapted to the IV curve to give Rp. ¹⁷

Further parameters that can be determined directly from the measurement data are the relative permittivities e of the individual organic layers. First, one has to ensure that the C-f plateau at V _DC = 0 V corresponds to the geometrical capacitance of all the semiconducting layers. This was checked by multiplying the capacitance with the total thickness. All devices with varied NPB and NBPhen layer thicknesses give almost the same value which indicates that we are probing the full device. Especially if highly doped or intrinsically polar layers are included, this condition is not necessarily fulfilled and the extraction of the relative permittivity is more difficult.19,21 Here, the geometrical capacitance of the total device is determined from the C-f plateau at V _DC = 0 V. Assuming a parallel plate capacitor for each layer, the total capacitance is given by a serial connection of each layer’s capacitance. Thanks to the variation of the NPB (NBPhen) thickness we can directly plot C-1 vs. d _NPB (dNBPhen), which is supposed to be linear with a slope corresponding to the permittivity of NPB (NBPhen). The permittivity of TCTA and the EML cannot be independently extracted as there is no thickness variation of these layers.

In order to perform a full electro-optical simulation we also need refractive indices of all the layers. Most of them are available from the Setfos database or were taken from literature. ^59,61,62 For the EML spectroscopic ellipsometry was performed on a film. Furthermore, the emitter source spectrum was measured by PL on a pure emitter film. The HOMO and LUMO energies were taken from literature. ^{16,20,22,63–65} Due to the uncertainty in the determination of the LUMO, ^66,67 this parameter was allowed to be changed by up to 0.5 eV during the fitting. The remaining electrical parameters, i.e. hole and electron injection barriers, density of chargeable sites, hole and electron mobilities for each of the layers were guessed or taken from literature if available. ^68–71 Here, the respective mobilities are described by the field-dependent Poole-Frenkel mobility model,

where μ ₀ is the zero-field mobility and γ the field-enhancement factor.By the choice of this mobility model, with only two free parameters, we follow the principle of Occam’s razor. This parametrization of charge mobility does not include any direct charge carrier dependence. If temperature-dependent measurements would be simulated, both charge- and field-dependence could be considered more easily by using the extended Gaussian/correlated disorder model (EGDM/ECDM). ^72,73 While the solution of the EGDM/ECDM has been presented before in steady-state, the implementation to solve multilayer structures in time and frequency domain is much more challenging. As the inclusion of transient and impedance data plays a key role in this analysis we choose the Poole-Frenkel mobility model which is the simplest model at room temperature that can well describe the presented experimental data.

Whenever more complicated device structures, doped layers or new transport materials are investigated, it is generally useful to first measure and simulate elementary devices, such as single-carrier devices, and then gradually increase the complexity of the investigated stack. ¹⁵ It has to be noted here, that the choice for the mobility values of the minority carriers, i.e. the electron (hole) mobilities in the NPB and TCTA (NBPhen), has a negligible effect on the simulated curves. Therefore, they are not used as fitting parameters in the following. To couple the electrical characteristics with the optical measurements, excitonic parameters have to be known. For our TADF system, radiative, non-radiative and intersystem crossing rates were extracted from a transient PL experiment with an emitter film sample. These parameters were used as input for the initial guess simulation. The full list of parameters used and obtain during our analysis can be found in our JAP publication (https://aip.scitation.org/doi/full/10.1063/1.5132599).

First, we notice that the simulated IV onset voltage (V_onset) is not in agreement with the measurements. The only parameter which influences this voltage is the built-in potential. In order to match the measured onset, we therefore decrease all LUMO levels by 0.15 eV in energy while keeping the energy barriers the same. After this adaption, the energy levels are now fixed for all the following optimization steps.

As mentioned above, thanks to the thickness dependence we already know that some hole mobilities have to be increased in the simulation to improve the agreement with the measurement in TEL rise and injection-CELIV. However, we cannot easily distinguish between the hole mobility in NPB, TCTA or in the EML. In order to quantitatively understand which parameter is influencing the electrical device characteristics in which way, we perform individual parameter sweeps starting from the initial  simulation shown in Figure 4. As examples we show the influence of the NPB, TCTA, EML zero-field hole mobilities, the field-enhancement coefficient for NPB as well as the zero-field electron mobilities in the EML and the NBPhen on the IV (Figure 5), the TEL rise, the injection-CELIV signal and the C-f plot (Figures S4 – S6). In order to simplify the graph, we only show the parameter variation for device S42.

The largest influence on the IV curve is observed for the zero-field hole mobility and the field-enhancement factor of the NPB. While the former gradually increases the forward current, the latter also influences the power law of the current above 2.5 V. Thus, we will clearly need to increase the latter in order to fit the experimental data. Interestingly, the charge carrier mobilities in the EML mostly affect the region just after V_onset. Especially the electron mobility should be decreased to match the slope of the IV between 2 – 3 V. The TCTA hole mobility and the NBPhen electron mobility have only a minor effect on the simulated IV curve. Using only these 6 parameters would allow to obtain a decent fit of the IV curve. However, we also want to include the impedance and transient signals in order to increase the reliability of the extracted parameters. This is illustrated by the cases of the hole and electron mobilities in the EML which affect the IV only slightly. While the EML hole mobility clearly shifts the TEL onset, the electron mobility has no effect on the onset time (Figure S6). In the case of injection-CELIV, both mobilities have an influence on the peak position t_max: increasing or decreasing the EML hole or electron mobility shifts the peak position to shorter or longer times, respectively.  This example highlights the importance of combining several measurement techniques, as the unknown parameters are reflected differently in each experiment.

As a next step, we will adjust the parameters to match the simulation with the selected key measurements. This can either be done with an automatic optimization/fitting algorithm or by manual tweaking of single parameters according to the observations in the parameter sweeps (Figure 5 and Figures S4 – S6). Automatic least-square fitting algorithms, such as Levenberg-Marquardt or simulated annealing, try to minimize the error χ between the measurement m and the simulation s for N user-defined target points which can be weighted by a factor w.

Due to the multiple devices and measurement techniques which should be optimized together, the selection of the target points and individual weighting factors w_i is not straightforward. The optimization is therefore done manually in this example. Furthermore, the fit quality is not evaluated quantitatively by an absolute error but judged by eye. The best fit that was obtained for the current model (no traps included) is shown in Figure 6.

This first fit is already able to reproduce the IV curve above turn-on very well. Also, the trends with varied NPB and NBPhen layer thickness in the TEL and injection-CELIV curves are well reproduced. However, the injection-CELIV peak is clearly too large which indicates that too many charge carriers are accumulated and extracted in the modeled device. The simulated impedance data shows quantitatively the same result as the measured C-f and G-f characteristics at 3 V. The NBPhen thickness dependence is however not correctly reproduced. Finally, also the fit of TEL decay (Figure 6(f)) still needs some fine tuning of the parameters. This measurement is the only one which is dominated by excitonic parameters, i.e. K_rad,s, K_nr,T, k_isc and K_risc. These parameters can therefore be independently adapted to match the TEL decay signal. The resulting simulated TEL decay curve is shown in Figure 8(f).

The other electro-optical simulations cannot be adapted further to improve the agreement with the experimental data using the existing model that does not consider any trap states. Therefore, according to the workflow in Figure 2, we shall diversify the model and include additional material parameters. The main features that we want to reproduce in a better way with the simulations are (1) the injection-CELIV peak area which is currently too large in the simulation, as well as (2) the NBPhen thickness dependence in the C-f and G-f data. As mentioned before, the injection-CELIV transient at V_offset = 3 V is mainly determined by hole transport parameters. We therefore first analyze the effect of including hole traps in the NPB, TCTA and the EML layer. We find that including hole traps in the EML does not affect the injection-CELIV peak at all. Therefore we do not further consider any hole traps in this layer. In striking contrast, for the NPB and TCTA layers, the hole trap density as well as the trap energy in TCTA have a significant impact on the CELIV peak area. The corresponding simulations of the injection-CELIV and the IV curve are shown in Figure 7. In more detail, we see that both the TCTA hole trap density and the trap depth have a very similar influence on the CELIV transient and to some extent also on the IV curve. This indicates that the parameters “trap density pt” and “trap depth ΔE_t,p” are correlated, which will also be discussed and quantified in Section III. D. Furthermore, we see from Figure 7 that including as many TCTA traps as necessary to match the injection-CELIV peak area will clearly deteriorate the IV behavior. We will therefore also need to include some traps in the NPB layer.

The effect of the NPB trap density is shown in Figure 7(c,f). As for the TCTA layer, the traps’ depth has a very similar effect on the simulated curves as their density (see Figure S7). In contrast to the traps in the TCTA layer, which mainly shift the CELIV peak in time while the qualitative shape is maintained, the traps in the NPB  also change the shape of the injection-CELIV peak. In some cases, even a double peak is observed. In order to find a good fit for both the IV and the injection-CELIV data we need to include and optimize the trap parameters in both the NPB and the TCTA layer. As an additional positive effect, the inclusion of hole traps in TCTA and NPB leads to a small TEL overshoot in the simulation which was also confirmed experimentally. Without trap states, the simulated TEL signal never exhibits an overshoot. The overshoot height can be adjusted by the capture rate of hole traps in the TCTA and NPB layer.

Before the reliability of the extracted material parameters is analyzed in the next section, we show how the simulation can further be used to understand device performance and analyze possible loss mechanisms. First, it is possible to analyze electrical profiles of the device at any working point. As an example, Figure 9 shows the recombination profile of the S42 device as a function of voltage. The recombination between 2 – 3 V is almost constant in the full EML. With higher voltages the recombination zone gets narrower and settles at the EML/NBPhen interface. This is in agreement with previous findings by Peng et al., who reported an exponentially decaying emission zone at the ETL interface with a width of about 7 nm using a very similar OLED structure.²² Such an accumulation of charge carriers at one interface of the EML is not ideal because degradation of organic molecules will predominantly occur at locations of high charge carrier and exciton concentrations.^22,74,75 A broadening of the emission zone is therefore desirable to increase the stability of this OLED structure. In this example, the electron transport in the EML host-guest layer is assumed to occur on the emitter molecules themselves. Because the emitter concentration is low (20%) the electron transport is limited, which is also reflected in the extracted mobility parameters for the EML. The zero-field electron mobility is about seven orders of magnitude smaller than the zero-field hole mobility of the EML. Balancing the EML mobilities would broaden the recombination zone which could improve the lifetime of these devices. Such an electron mobility increase can, for example, be achieved by using a co-host approach or by increasing the emitter concentration.²²

A second criterion to assess the reliability of the extracted parameters is the sensitivity of the investigated characterization techniques to the material parameters. Parameter sensitivity analysis can be divided into a local and a global one. The first evaluates the influence of a parameter variation around a specific local working point on the agreement between simulation and measurement. Here, the sensitivity is evaluated for the final fitting parameters as working point, yielding the sensitivity of the final fit (Figure 8) on each of the parameters.

Global sensitivity analyzes the combined effect of varying all parameters within a certain range. These randomized simulation results are typically shown as scatter plots. ⁷⁶ While it can give more general insights into the sensitivity of specific measurement techniques to material parameters, it is more difficult to be quantified compared to the local sensitivity analysis and computationally more demanding. In the scope of this tutorial we therefore focus on the local parameter sensitivity at the final fit.

In previous publications on single-carrier devices and on organic as well as perovskite solar cells we have studied the qualitative influence of the involved parameters by a systematic variation around the values from the final fit. ^19,27,60 In the present paper on OLEDs we show the sensitivity of the electron mobility parameters µ_0,e and γ_e of the EML as well as the zero-field hole mobility in the NBPhen layer on the selected measurements (Figure S9). Due to the large number of parameters, devices and performed measurements this will not be presented here in detail. Both EML parameters show a certain influence on one or the other experiment. Interestingly, they seem not to influence the TEL rise (Figure S9 (e, k)) while both have for example an effect on the IV curve between 2 and 3 V (Figure S9 (a, f)). Therefore these parameters cannot be accurately determined from a TEL rise fit alone. This finding again demonstrates how important it is to combine different measurement to be sensitive to as many parameters as possible. In Figure S9 we also show a parameter sweep of the NBPhen hole mobility. Note that this parameter was not used as a fitting parameter. A variation of ± one order of magnitude does not have any influence on the selected measurement techniques. One can say that we are not sensitive to this parameter and one should thus not try to determine this parameter from the selected experiments.

In order to obtain a more quantitative analysis of the sensitivity with respect to the involved model parameters, we again consider the error f_i between the measured data and the simulation at N different user-defined target points (see equation (5)) and define the local sensitivity as

The calculated local sensitivity can take values between zero and infinity. The results are shown in Figure 11 and ordered from high to low. There is no general definition of a value that is still acceptable. In order to simplify the intepretation of the calculated local sensitivity, three groups of “high”, “medium” and “low” sensitivity were defined according to the observations in the parameter sweeps (Figure S9). In the following discussion we will show three examples (marked in blue in Figure 11) to illustrate the different sensitivities.

The hole (minority carrier) mobility in NBPhen is an example of “low” sensitivity. Such parameters cannot be determined accurately from the global fit. In this study, the minority carrier mobilities were not used as a fitting parameter. The field-enhancement coefficient of the EML electron mobility is a good example for “medium” sensitivity. It has some effect on some of the selected experiment but its influence is not that large. Parameters with a “high” sensitivity show a clear trend in the experiment when they are varied. Those parameters can safely be determined from the global fit, provided that they also show a low correlation to other parameters. Overall, we find that owing to the combination of 4 different experimental techniques and 8 different devices, the model is sensitive enough to all fitted parameters.