Universal approach for diffusion quantification applied to lead halide perovskite single crystals

A universal approach to calculating diffusion coefficients in lead halide perovskite single crystals, which have ionic and mixed ionic–electronic conductivity, is proposed. Using impedance spectroscopy, it is demonstrated how to model a non-ideal Warburg element and transmission line equivalent circuit to identify ionic diffusion in the material. The proposed method is applicable to samples of any thickness and electrical properties. Additionally, it is shown how to overcome the challenges of low-frequency impedance measurement and the non-ideal behavior of the elements through extrapolative modeling and approximation.


Introduction
Understanding charge carrier dynamics in the material can be crucial for deep comprehension and handling of degradation and operational processes of electronic and optoelectronic devices. Lead halide perovskites, apart from being a semiconductor, additionally show pronounced ionic conductivity [1,2], which is the result of ionic diffusion, and thus a non-uniform concentration of charged particles in the material. This type of charge transport appeared to be responsible for many unique properties of lead halide perovskite that had been long debated and various explanations had been suggested [3][4][5][6]. For instance, in [7,8], it was proposed that capacitive anomalies at low and intermediate frequencies can be associated with ionic migration in methylammonium lead bromide single crystals.
Widely reported and discussed [8], hysteresis in perovskite solar cells, despite many various suggestions in the past, [9] at the present moment is mainly associated with slow ionic charge transport [10,11]. Several studies have reported [13][14][15][16] that halide perovskites containing IV-group B-cations have a discernible ionic contribution to the overall electrical transport characteristics.
Impedance spectroscopy investigation on lead halide perovskite films [17][18][19] revealed the presence of a low-frequency contribution, which was attributed to the migration of slow-moving charged species. Moreover, Warburg-like semi-infinite diffusion at low frequencies was observed in lead halide solar cell-like architecture under light illumination. This behavior is typical for mixed conductors [20], which is indirect evidence of moving species in the material. The ion migration hypothesis was supported by a large amount of evidence from detailed research. Ionic defect conductivity is well known in similar materials like metal oxide perovskites. Moreover, it was discovered that a cubic or orthorhombic shape correlates with high ionic conductivity [21,22].
Following the aforesaid, charge carrier dynamics in single crystals is of utmost importance, as it corresponds to the intrinsic property of the material that helps to understand various processes and operational mechanisms [23,24], such as operational and/or temporal stability [25,26] of optical and electrical devices. In this regard, several great efforts have been done to explain charge carrier dynamics in perovskites [27,28]; however, some of them were performed using complicated layered structures that can distort diffusion coefficient due to the interfacial influence, while others have limited applicability due to the complicated sample preparation. Moreover, in the case of the measurements provided on thin films, given an extremely long ionic diffusion range, the reduced dimensionality of samples dramatically affects the charge carrier dynamics in the material, requiring careful consideration in quantifying diffusion [29][30][31].
Here, we present a simple and universal approach for diffusion quantification applied to lead halide perovskite single crystals. Single crystals, prepared by inverse temperature crystallization [32], were used as they show exceptionally low intrinsic defect density [33] and generally are easy to synthesize. Nevertheless, we assume that the strategy described below allows direct monitoring and therefore quantifying of the ionic-related conductivity regardless of the presence of electronic current, and is applicable for an arbitrarily mixed conductor possessing ionic diffusion.
The material's resistive, capacitive, and inductive properties all contribute to its electrical response in the frequency-time domain. By understanding these response characteristics, we can explain various processes that are happening within the material. For this reason, impedance spectroscopy is a powerful tool that measures the electrical response, depending on the material's electrodynamic properties [34,35]. By monitoring complex impedance values and phase shifts at various frequencies, this technique allows separating resistive, capacitive, and inductive current contributions by using appropriate mathematical modeling.

Warburg-like diffusion model
While the dynamics of ionic processes are much slower than electronic ones, their characteristic features tend to appear at relatively low frequencies. Therefore, by locating characteristic turnover frequency, it is possible to quantify diffusion in the material. The model for diffusion quantification described by Peng et al. [28] fits well, for thin films when the sample thickness is less than the diffusion length. Lead halide perovskites are well known for their exceptionally long diffusion length in the range of micrometers [36]; thus, the above-mentioned model can be successfully applied for the films used for most optoelectronic devices, as their thickness rarely exceeds 0.5 um. However, for thicker films or single crystals, the model is not valid once the thickness of the sample is larger than the diffusion length. Moreover, this technique is barely applicable to measure diffusion coefficients in materials with shorter diffusion lengths, e.g.. organic semiconductors [37].
In this paper, a sandwich structure was used, and the Warburg model of diffusion for flat parallel electrodes can be applied. A Warburg impedance element Z W is nearly always associated with a charge-transfer resistance and a double-layer capacitance [38]. Therefore, once Warburg-like characteristics of the impedance response are observed, the diffusion coefficient can be determined. However, it can be difficult to identify a Warburg element in practice due to two main issues: non-ideality and the low-frequency range. In this regard, we denote non-ideal diffusion as 'Warburg-like diffusion'; later, we discuss why such diffusion is problematic to detect and evaluate.

Diffusion non-ideality
Mathematically, the Warburg element total impedance at equilibrium can be defined by the equation (1) [39], as the magnitude of the Warburg impedance is reversibly proportional to the square root of angular frequency according to: where Warburg coefficient is related to the concentrations and diffusion coefficients of the oxidized and reduced species ( C ox , C red , D ox , and D red , respectively) by: where R is the gas constant (J/K/mol), T is the temperature (K), n is the number of exchanged electrons, F is the Faraday constant ( C/mol ), A is the area of the electrode ( m 2 ∕s ), C ox and C red are the bulk concentrations of electroactive species in their oxidized, resp. reduced form ( mol/m 3 ), and D ox , and D red are diffusion coefficients of oxidized, resp. reduced species ( m 2 ∕s ).
Thus, the Warburg circuit element is typically associated with the existence of a linear relationship with a slope of -1/2 on the log of a Bode plot ( log |Z| versus log ( ) ), and it exhibits a 45 • phase shift. Moreover, in the case of ideal Warburg behavior, at the frequency range where the phase shift is frequency independent, plotting both real (R) and imaginary (X) part of the impedance Z versus 1∕ 1∕2 results in two straight parallel lines where the slope is directly the Warburg coefficient . Warburg pseudocapacitance ( C W ) cannot be regarded as a pure capacitor, because it is dependent on faradaic charge transfer across the interphase.
However, in the original study [39], Warburg reports smaller phase angles for the different electrodes used in the experiment; unfortunately, at the same time, the frequency range is not reported. A phase angle of 45 • is a result of deviation from a semi-infinite diffusion condition. Nevertheless, one can notice that experiments on real-life samples or devices (solar cells, batteries, detectors, etc.) often show non-ideal Warburg-like diffusion with phase angles, more often higher, or sometimes lower than -45 • . Lower phase angles are associated with, for example, additional leakage currents that make the whole system less capacitive, while higher angles are a result of factors that add additional capacitive behavior, such as extra porosity, which can charge, like a capacitor, but will not be a part of the diffusion component.
In cases where diffusion is expected, such as lead halide perovskites as discussed in this paper, the diffusion model cannot be considered an ideal Warburg model. As a result, the dependencies of R and X on 1∕ 1∕2 will not be linear or parallel, making it difficult to accurately estimate the diffusion coefficient.
As a result, the presence of diffusion may go undetected because the usual methods for identifying Warburg-like behavior will not be applicable. To detect diffusion in these cases, we propose the following approach: • Find a low-frequency interval, where the phase shift is between 0 • and -90 • and frequency independent. • Given the frequency interval meeting the previous condition, find a slope of ln(Z) vs ln( ) , which gives a value of a so-called ideality factor k. • Plot both R and X vs. 1∕ k , which we denote as a modified Warburg impedance plot.
In this case, if the plots in the mentioned frequency range are linear, diffusion can be quantified. Note that R and X may not be parallel, which will be discussed later. However, it is worth noting that R and X may not be parallel, which will be discussed further on. As a result, the slopes for the real and imaginary parts of the impedance, R and X , respectively, are distinct. Thus, total non-ideal Warburg-like impedance at the equilibrium can be defined by the equation (3): where R and X are non-ideal Warburg coefficients for resistance (R, Z R ) and reactance (X, Z X ), respectively. From equation 3, both real and imaginary parts of the impedance can be derived as follows: Thus, by plotting ln (Z R ) or ln (Z X ) vs ln ( ) (i.e., modified Bode plot for the real and imaginary parts of impedance as ln Z X = ln X + k ln( ). a function of angular frequency), the ideality factor k can be easily determined from the slope of this dependence; moreover, from the intersection on the y-axis, the Warburg coefficient (representing the ideal case of the Warburg element) can be obtained at the same time.

Diffusion detection frequency
At a constant temperature, the diffusion coefficient mainly depends on the intrinsic characteristics of the material that define the diffusion coefficient [39]. The relation between the diffusion coefficient and the turnout frequency range on the impedance spectrum was previously discussed by Bisquert [40] and Janssen et. al [41]. Therefore, at the impedance spectrum, diffusion may occur in an extremely lowfrequency range, micro-Hz or even nano-Hz. In an ideal case, the frequency interval where diffusion occurs should be in a range where the phase shift is negative and frequency independent, -45 • as can be derived from equation 1. In nonideal diffusion conditions, higher phase angles between 0 • and -45 • (more resistive behavior) may be caused by contact/interface non-idealities leading to charge carrier dispersion, while lower phase angles between -45 • and -90 • (more capacitive behavior) may occur in porous materials where pores act like capacitors that accumulate charges, but are not part of the diffusion component.

Warburg diffusion element equivalent circuit
As discussed above, using impedance spectroscopy for diffusion detection/quantification can be quite problematic, as the diffusion can appear at an extremely low-frequency range; therefore, measurements would last weeks and even longer, which is practically often not acceptable. However, the power of impedance spectroscopy lies in the ability to extrapolate using equivalent circuit modeling [35,42]. Fortunately, equivalent circuit for both Warburg diffusion [43][44][45][46] and lead halide perovskites [18,[47][48][49] are well reported and described in the literature. Therefore, it is possible to get a hint on the equivalent circuit of the whole system. In the present paper, to fit impedance response within the frequency domain, an equivalent circuit, depicted in Fig. 1, was used. The diffusion element was modeled according to the literature [50][51][52][53], and the complex impedance Z dif f of the electrical circuit including the Warburg diffusion element can be described by equation (6): where is the Warburg coefficient, Z 0 is the magnitude of the impedance at = 1 rad/s , is the angular frequency, j is the imaginary unit, and R 3 and C 3 are the low-frequency resistance and capacitance, respectively, while the full impedance of the whole circuit can be defined as shown in equation (7): where R 1 and C 1 are high-frequency resistance and capacitance, R 2 and C 2 are middle-frequency resistance and capacitance, Z dif f is the low-frequency complex impedance of the circuit consisting of the Warburg element, R S is the series resistance, and L S is the series inductance. Therefore, using the above-described relations for the full impedance Z, one can fit the impedance response in the frequency domain.
Fitting precision for the current measurements is shown in supporting materials.

How do we quantify diffusion?
Given the above, we propose a protocol for Warburg-like diffusion quantification: • Ideality factor k can be extracted correctly only if we know the behavior at very low frequencies, where the diffusion takes place (often nano Hz). To obtain k we need to fit the model to the existing data, using the equivalent circuit given in Fig. 1 defined by the correspondent equa-tion6, resp. 7. • Once we fit the model on the measured data, we are able to extrapolate the data down to nano Hz, to detect the potential diffusion frequency range where the phase angle is negative and frequency independent. • From the slope of ln (Z) vs ln (f ) or ln ( ) at the frequency range obtained from above, the k coefficient can be obtained. It has to be noted that k can be obtained also from the phase angle, at the diffusion frequency range. • The ideal case assumes linear regression on the R − X plot vs. −−1∕2 , where both resistance and reactance dependences are straight and parallel with respect to −−1∕2 . However, in the case of non-ideal Warburg behavior, R and X are not parallel and non-linear, so the coefficient will depend on the frequency, and there will be a significant difference between R and X resulting in the interpretation ambiguity once again. • However, if we plot R − X vs. −−k , it linearizes both dependencies and makes R and X almost equal.
After deriving R and X values, the diffusion can be easily quantified.

Results and discussion
In the ideal case, Warburg diffusion can be modeled using the so-called 'ideality factor' k = 1∕2 . Nevertheless, in 'real-world' samples having such a value for k is rather an exception. Usually, one can observe lower values because of any sort of leakage currents, and also higher values can be observed in the case of porous materials. In these cases, the Warburg coefficient is problematic to determine, as nonideality results in a non-linear R − X plot vs. 1∕ 0.5 . We propose that in such cases, ideal factor 0.5 should be replaced by k factor; therefore, the plot should be represented as R − X plot vs. 1∕ 0.5 . As an example, let us analyze the data obtained from the MAPbI 3 perovskite monocrystal sandwich structure. In the first step, the diffusion frequency range has to be located. As we discussed above, in this range, the phase angle should be negative and frequency independent (Fig. 2b), the slope of ln (Z) vs ln ( ) in this frequency range (Fig. 2a) gives us the ideality factor k = 0.21 , i.e,. non-ideal Warburg behavior, and phase shift corresponds to approx. 21 • at low frequencies.
After introducing k = 0.21, one can obtain both R and X modified Warburg constants. By the way, as described above, we can get the Warburg constant directly from the Bode plot according to equations (4) and (5). To quantify the Warburg coefficient, we finally used the mean value of R and X . As a demonstration, for MAPbI 3 perovskite single crystal with ideality factor equal to 0.21, the calculated Warburg coefficient (avg value) from the modified R − X plot vs. −−0.21 is approximately k = 2.29 × 10 9 ∕s 0.21 ( R = 3.43 × 10 9 ∕s 0.21 , X = 1.14 × 10 9 ∕s 0.21 ) . T h e determined Warburg coefficient (average value) from the modified Bode plot is approx. = 2.49 × 10 9 ∕s 0.5 ( R = 3.73 × 10 9 ∕s 0.5 , X = 1.24 × 10 9 ∕s 0.5 ).
To find a correlation between an ideal and a non-ideal Warburg behavior, we tested a large auto-generated set of models (for more details see Supporting Information) with the equivalent circuit described above in Fig. 1. A general trend between an ideal and a non-ideal Warburg component that follows a quadratic function, empirically defined by an equation y = 1.647x 2 − −1.568x + 1.373 , is shown in Fig. 3a. Due to this relationship, the ideal Warburg element (and thus the final Warburg coefficient) can be easily estimated even if we have non-ideal Warburg behavior in our system. However, what arises here is the Warburg coefficient unit, due to the different powers/roots for angular frequency in the graphs. In this regard, it is necessary to introduce a new correction factor corr = 0.5 ∕ k rad 0.5−k s k−0.5 , which adjusts the finite unit of the Warburg coefficient (especially for purposes of non-ideal Warburg behavior) so that it is correct for evaluating the diffusion coefficient with its unit m 2 ∕s.
As an example, according to these findings, the abovedetermined Warburg coefficient from the modified R − X plot ( k = 2.29 × 10 9 Ω∕s 0.21 for non-ideal Warburg behavior with the ideality factor equal to 0.21) has to be multiplied Fig. 2 A The modified Bode plot, points for real measured data, solid lines for modeled data up to 1 μHz, dashed lines for interpolation of data in low frequencies; B a phase shift between the real and imaginary part of Z plotted against ln ( ) with a line indicating the area where the phase shift is independent of frequency; C a plot of the non-ideal Warburg-like impedance; D the modified plot of the nonideal Warburg-like impedance using coefficient k for ideal or nonideal behavior Fig. 3 A Correlation between an ideal and a non-ideal Warburg component. B The Warburg coefficients for a Warburg diffusion element equivalent circuit ( R 3 = 5GΩ , C 3 = 1nF , Z W = 1GΩ ) with variable ideality factors and its slope which follows the standard value for an ideal Warburg component by a correction factor of 1.12, and after that, the 'ideal' Warburg coefficient equals = 2.57 × 10 9 Ω∕s 0.5 . Here, we finally used the mean value of the determined Warburg constants (modified Bode plot, modified R − X plot) to quantify the diffusion. The final obtained Warburg constant for MAPbI 3 perovskite single crystal with an ideality factor of 0.21 is thus equal to = 2.53 × 10 9 Ω∕s 0.5 . Once we are able to calculate the Warburg coefficient as described above, the diffusion coefficient can then be easily quantified according to the modified equation (2) as follows: Note that other parameters from the equation, such as the dielectric constant of the material used and doping (defect) density calculated from Mott-Schottky analysis, were determined by impedance spectroscopy [28]. Final parameters obtained for all the measured crystals of different types ( Cl − , Br − , I − ) are summarized in Table S1, Table S2, and  Table S3 (see Supplementary Material). Consequently, we determined the diffusion coefficients for all three types of macroscopic single crystals: D = (2.68 ± 1.94) × 10 −8 cm 2 ∕s for MAPbCl 3 , (1.12 ± 0.61) × 10 −6 cm 2 ∕s for MAPbBr 3 and (3.63 ± 2.36) × 10 −9 cm 2 ∕s for MAPbl 3 . Comparing these results with the ones published in the literature, we have to mention that the reported diffusion coefficient estimation differs vastly depending upon the measurement technique and sample preparation [54][55][56][57][58][59][60]. The overall diffusion coefficients for MAPbI 3 are comparable with the ones obtained by transient ion-drift measurements in thin films [55][56][57]. At the same time, both much higher and much lower values can be obtained depending on the method of measurements and sample preparation procedure. Transient absorption experiments using ultrafast microscopy [61] usually result in higher values for a diffusion coefficient, as it takes into account a contribution from both electrons and holes, while photoluminescence-based [27,[62][63][64] and electrical [28] measurements usually underestimate diffusion coefficient and diffusion length. Moreover, in the case of electrical measurements of lead halide perovskites, the presence of moisture and ambient conditions play a major role and, therefore, measurements obtained using non-encapsulated samples can differ radically.

Conclusion
In this study, we present a universal approach to quantify diffusion. The method is applied to lead halide perovskites, a material with mixed ionic-electronic conductivity that is widely used in electronics and optics. Although diffusion quantification using impedance spectroscopy has been previously described, there are several challenges that can arise in obtaining well-defined diffusion coefficients for the material, such as a material thickness larger than the diffusion length, non-idealities present in real-world systems, or a low-frequency range of the diffusion that makes it difficult to accurately locate the turnover frequency on the Nyquist plot. We describe how to overcome these challenges and accurately quantify the diffusion coefficient. It is worth noting that, while the described approach is specific to lead halide perovskites, it is also applicable to any other type of material.

Experimental
In this study, we applied our universal approach to three types of lead halide perovskite monocrystals: MAPbCl3 , MAPbBr3 , and MAPbI 3 synthesized by inverse temperature crystallization [32]. All three types of crystals, with various thicknesses, were analyzed to obtain reliable results and demonstrate the reproducibility of the method. The thickness of the crystals (the distance between the two electrodes) was in the range of millimeters, which was greater than the diffusion lengths in all cases. Additional information about the crystal synthesis and measurement equipment can be found in the supplementary information section.

Acknowledgements
The authors thank GACR project 20-29499Y by Czech Science Foundation, for the financial support.
Funding Open access publishing supported by the National Technical Library in Prague.

Conflict of interest
The authors declare no conflicts of interest or personal relationships related to the work reported in this paper.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.