Abstract

Potentialities of density functional theory (DFT) based methodologies are explored for photovoltaic materials through the modeling of the structural and optoelectronic properties of semiconductor hybrid organic-inorganic perovskites and GaAs/GaP heterostructures. They show how the properties of these bulk materials, as well as atomistic relaxations, interfaces, and electronic band-lineups in small heterostructures, can be thoroughly investigated. Some limitations of available standard DFT codes are discussed. Recent improvements able to treat many-body effects or based on density-functional perturbation theory are also reviewed in the context of issues relevant to photovoltaic technologies.

1. Introduction

Photovoltaic (PV) solar electricity is one of the key technologies of the 21st century to reduce the world’s reliance on fossil fuels for energy generation. Reduced costs and higher conversion efficiencies are of crucial importance to make PV-based technologies economically more competitive. The quest for quality and improved performances for future solar cells has attracted a tremendous research effort over the last decade towards the semiconductor heterostructures, nanostructured materials, and thin films. Various approaches are used ranging from high-cost/high-performance III-V technology, multiple junctions, and concentrator systems to the low-cost thin films [1]. Obviously, the design of novel and/or efficient PV devices requires a deep understanding of underlying material’s properties including chemical composition, mechanical, electrical, and optical properties, which could be achieved from state-of-the-art ab initio approaches. Such knowledge is also desirable to reach low-cost PV cells composed of earth-abundant elements based materials capable of full recycling.

Among the large panel of available theoretical approaches, the density-functional theory (DFT) has become overwhelmingly popular. This success greatly relies on the fact that no input adjustable parameters are needed and that efficient numerical codes exist and a strong versatility especially in the description of the ground state properties of semiconductors and metals. Increase in computing power has afforded further capabilities in system’s size that DFT methods can handle. However, the current limit does not yet reach the window 10000–10 million atoms involved in the active device region of PV cells. Indeed, modern semiconductor optoelectronic devices have a feature size of few nanometers, including nanostructures like quantum wells [2] or quantum dots [3]. Such systems consist of involved two- (2D) and/or three- (3D) dimensional geometries composed of multiple materials and/or alloys. Nevertheless, DFT methods are useful to get insight into physical phenomena of each component part of the device, that is, separate materials or small heterostructures. For example, a quantitative design of PV requires a reliable prediction of the electronic band-gaps, band-lineups, and effective masses.

It is well known that DFT based on the local density approximation (LDA) [4] or generalized gradient approximation (GGA) [5] does not reproduce accurately the excited states of compounds. Conversely, hybrid methods that include a fraction of Hartree-Fock exchange may circumvent the band-gap problem, but their performances strongly depend on the material of interest. Interestingly, the HSE06 hybrid functional proposed by Heyd et al. [6] is a good alternative to compute band-gaps, band offsets, or alloys properties [7, 8]; even so, it fails to reproduce the direct-indirect crossover in GaAsP alloys [8]. Alternatively, reliable results can be obtained from many-body perturbation theory (MBPT), especially within the GW approach (GW, where G stands for Green’s function and W for the screened potential) that may be used in a perturbative scheme [911] or self-consistently [12, 13]. For charge separation, key quantities are the ionization potential and electron affinity, respectively, of the donor and acceptor materials, as they control the relative alignment of electron and hole levels. The drop of the interface’s potential in heterostructures can be estimated by DFT within a superlattice (SL) approximation, which consequently offers a nice estimate of the band-lineup. For a complex stack, a GW treatment is beyond reach, but the DFT potential drop at the interface can be efficiently corrected by the GW eigenvalues obtained for the bulk valence band states [11, 14, 15]. A good estimation of alloys electronic properties is also a difficult task in DFT simulation. Indeed, semiconductors, even conventional ones, undergo a significant bowing of the band-gap; that is, the band-gap energy follows , where is a bowing parameter. Supercell approaches for statistically random alloys are more adequate than virtual crystal approximations. For special alloy compositions, equivalent results can be obtained from DFT with special quasirandom structures (SQS) (small supercells) that reproduce mixing enthalpies and atomic correlations of very large supercells [8]. Noteworthily, SQS models include chemical mixing, strain, and atomic relaxation effects.

In addition, semiempirical methods for mechanical or electronic properties studies, as the valence force field (VFF) and the tight binding approximation [16] or elasticity and the k.p method [17], require accurate electronic parameters as input. They can be derived from DFT calculations or experiment. For example, quantitative estimates of electromechanical tensors of bulk materials are achieved from density-functional perturbation theory (DFPT) [1820]. An efficient use of the “” theorem, where only by-products of a first-order perturbation calculation are required, provides second- and third-order derivatives of the total energy, provided that atomic-displacement variables have been eliminated. For second-order derivatives, various physical responses of insulating crystals can be obtained, including elastic constants, linear piezoelectric tensors, and linear dielectric susceptibility, as well as tensor properties related to internal atomic displacements like Born charges and phonons [20, 21]. For third-order components, related to physical properties such as nonlinear electrical susceptibilities, nonlinear elasticity, or photoelastic and electrostrictive effects, finite difference techniques and symmetry analysis must be associated with the DFPT method [22]. Finally, excitonic and transport properties can also be studied by DFT methods, but they are beyond the scope of the present paper.

This paper aims at illustrating, with two examples, some of the DFT’s potentialities for PV technologies. First, we will focus on hybrid organic/inorganic perovskites (HOP) that open new routes for solar cells (Section 3.1). In particular, a direct optical transition and isotropic activity are predicted for a model 3D system. Secondly, the investigation of structural, elastic, and electronic properties will be illustrated on III-V semiconductor heterostructures for GaP/Si pseudosubstrates (Section 3.2). Computational details are given in Section 2.

2. Computational Details

Total energy DFT calculations were carried out using the ABINIT code [20] within LDA or GGA. A plane-wave basis set with an energy cut-off of 340 and 680 eV for compounds investigated, respectively, in Sections 3.1 and 3.2, was used to expand the electronic wave-functions. The reciprocal space integration was performed over Monkhorst-Pack grids [23]. The energy was computed from the linear response method and convergence is accurately reached with tolerance on the residual potential which stems from differences between the input and output potentials. Pseudopotentials were constructed either from the Fritz Haber Institute (FHI) format [24] or using the Hartwigsen-Goedecker-Hutter (HGH) scheme [25]. The dielectric properties were studied from the density-functional perturbation theory (DFPT) based on the linear response theory implemented in SIESTA code [26]. All many-body GW calculations were performed at the G0W0 level [9]. The quasiparticles energies were converged with respect to the energy cut-offs, the number of k-points in the Brillouin Zone (BZ), and the number of bands used to compute the dielectric function and the self-energy. The plasmon-pole model was used to describe the dynamic dependence of the screening function [20].

3. Results

3.1. Optoelectronic Properties of Hybrid Organic/Inorganic Perovskites for High Efficiency and Low-Cost Solar Cells

Hybrid organic/inorganic materials have attracted increasing interest over the past decade due to their potential applications. In particular, hybrid nanomaterials are expected to offer a way towards enhanced performances of optoelectronic devices [27]. The versatility of the organic part affords the possibility of fine tuning material’s optoelectronic properties [28]. For example, it has been shown that the optical spectra of HOP can be easily tailored by varying the organic cation, which improves optical efficiencies and the tuning of the emission wavelength. Among them, self-assembled layered (2D) HOP have recently shown enhanced nonlinear optical properties in microcavities [29]. Moreover, following pioneering recent works [3032], 3D HOP based on relatively small organic cations have also been shown to drastically improve the photonic conversion in dye sensitized solar cells (DSSC) [3338] and expected to afford efficiencies up to 20% in the near future [39]. It is even predicted to open “a new era and a new avenue of research and development for low-cost solar cells … likely to push the absolute power conversion efficiency toward that of CIGS (20%) and then toward and beyond that of crystalline silicon (25%)” [40]. Our theoretical work on this rapidly evolving topic shows that these compounds can be considered now as a new class of semiconductors [41]. Indeed, compared to alternative strategies based on inorganic semiconductor quantum dots (QDs) or extremely thin absorbers coated upon the internal surface of a mesoporous TiO2 electrode, 3D HOP offer different benefits. Among those, one can note the ease and low temperature of synthesis, the tailoring of electronic properties by chemical substitution, hole transport, their high stability in dry air, and the formation of pn-like heterojunctions on TiO2.

DFT calculations on prototypes 3D or 2D HOP combined with symmetry analysis of the band edge Bloch states show that the ordering of the band edge states is found reversed compared to tetrahedrally bonded semiconductor structures [28, 41]. Moreover, from the computation of Kane’s energy parameters, we have explained the underlying mechanism of the transverse electric (TE) optical activity of 2D HOP resulting from a subtle interplay between the electronic structure and exciton binding energies [28]. Interestingly, it has also been shown that the spin-orbit coupling (SOC) induces a large splitting of the conduction bands, in comparison with those of the valence bands of cubic semiconductors. The importance of SOC was also put forward for 3D HOP [41]. To illustrate this point and the understanding gained from DFT methodologies for hybrid compounds, we consider one of the 3D HOP proposed for PV devices. It belongs to the CH3NH3PbX3 family (where X is a halogen atom) and exhibits a disordered average cubic phase (space group Pm3m) at room temperature. Structural disorder is associated with both the rotation of ammonium cations and tilt of the lead halide octaedra. CH3NH3PbI3 presents an ordered structure at low temperature, which is orthorhombic (space group Pnma), with a cell doubling when compared to the room temperature phase [41, 45] (Figure 1). The electronic degeneracy lifting is associated with the low temperature symmetry breaking and strain into the unit cell:  Ǻ,  Ǻ, and  Ǻ [45]. The low temperature cell strain is larger along the -axis. The theoretical study of 3D HOP predicts an isotropic absorption at high temperature for these materials used in the active zone of solar cells [41]. The triply degenerate CB without SOC of 3D HOP states is indeed associated with the vectorial representation of a cubic simple group. For the low temperature orthorhombic Pnma phase of the 3D CH3NH3PbI3 HOP, the number of CB states is doubled (Figure 2(a)). The fundamental electronic transitions for the CH3NH3PbI3 compound with and without SOC are found at the DFT-LDA level equal to 0.4 (Figure 2(b)) and 1.4 eV (Figure 2(a)), respectively. These values are a little bit smaller than the values computed at the DFT-GGA level [41]. Even though the fundamental transitions calculated without SOC for the low temperature Pnma structure (1.4 eV) compare nicely with the values obtained experimentally (1.5 eV, [30]), the agreement is fortuitous and stems from large error cancellations. SOC effect indeed strongly reduces the band-gap (Figure 2), without modifying the main character of the optical transitions [41]. At the same time, it is known that a DFT ground state computation systematically underestimates the excited states. This deficiency can be solved with the inclusion of many-body effects, namely, the GW self-energy correction for the band-gap and the subsequent resolution of Bethe Salpeter’s equation to account for the excitonic features. However, such calculations are far beyond available computational resources for large systems. In the corresponding double group including spinors, the CB of the high temperature cubic phase is split by SOC into twofold degenerate states and fourfold degenerate states [41]. A similar effect is predicted for the low-temperature phase (Figure 2(b)). The SOC splitting (on the order of 1.0 eV) is much larger than the one usually encountered in the valence band (VB) of cubic conventional semiconductors. The CB minimum is associated with the twofold degenerate and odd spin-orbit split-off (SO) states, leading for symmetry reasons to a strong and isotropic optical transition with the even VB states [41].

These findings differ from those derived for analogue 2D HOP with larger organic cations. In fact, the 2D densities of states obtained from 4F-PEPI ((C5H11NH3)2PbI4) [28] and that obtained from (C18H37NH3)2PbI4, whose structure is shown in Figure 3 [46], show quite different features. The band-gap remains located at the Γ-point are direct and associated with only three active Bloch states: a nondegenerate level for the valence-band maximum (VBM) and two nearly doubly degenerate levels for the conduction-band minimum (CBM1-2) ( Figure 4, [28]). Moreover, due to the 2D character of these layered compounds, no energy dispersion occurs along the Γ-X direction (characterizing the stacking axis in real space), illustrating the Pb-I-Pb bond breaking and the dielectric mismatch between inorganic and organic sheets. The fundamental transition of hybrid organic/inorganic layered (2D) perovskites is indeed expected to display a TE character, and it could be even further enhanced by orientational disorder introduced by the organic layer [28]. Contrarily, in the CH3NH3PbI3 HOP, an almost isotropic optical activity is expected from symmetry that should also be enhanced at room temperature with the disorder observed for the cubic phase [41]. The computation of the relative imaginary dielectric permittivity is computationally involved with such large systems. It is then limited in this work to a calculation without SOC, for the 3D HOP and for a 2D HOP compound with a shorter alkyl chain (C5H11NH3)2PbI4 [41, 47]. Relative imaginary dielectric permittivity computed at the GGA-level presented in Figure 5 for CH3NH3PbI3 (a) and (C5H11NH3)2PbI4 (b) confirms our earlier theoretical predictions based on symmetry analyses [28, 41]. The results are given for light polarizations along the three crystallographic axes (1–3) or for a polycrystalline sample in the case of CH3NH3PbI3 (a) and TE or TM polarizations for (C5H11NH3)2PbI4 (b).

3.2. III-V Semiconductors Heterostructures for GaP/Si Pseudosubstrates

High efficiency PV technology is based on III-V materials growth on expensive Ge or GaAs substrates [1, 48, 49]. With the combination of bulk GaInP and GaAs, 30% conversion efficiencies are attained. Higher conversion efficiencies than 40% can be reached using triple GaInP/Ga(In)As/Ge junctions in solar cells. Efficient tandem cells based on Si substrates have been recently proposed [5052]. They are based on alloyed heterostructures GaAsP(N)/GaP(N)/Si double-pn-junctions. The key point concerns a coherent growth of III-V semiconductors on a Si substrate. An initial deposition of thin GaP layer (with a lattice mismatch equal to 0.36%) and a subsequent overgrowth of the diluted-nitride GaAsPN compound lattice-matched to Si are a possible solution to obtain a 1.7 eV/1.1 eV tandem cell. Elaborated growth strategies have however to be used in order to overcome the problems related to extended structural defects formation like antiphase domains and MicroTwins [5052].

GaAsP(N)/GaP(N) material combinations are good candidates to reach optimal absorption at 1.7 eV [53, 54]. Accurate band structure calculations of strained Ga(NAsP) heterostructures have been performed within the framework of the extended-basis sp3d5 tight-binding model [16]. DFT calculations are helpful to yield comparisons and band parameters for semiempirical approaches. They are however limited for the simulation of ultrashort superlattices or binary alloys. Figure 6 shows the atomic positions of GaAs(4)/GaP(4) SLs grown along the direction. The GaAs layer is strained onto a GaP substrate where the DFT lattice parameter is of  Ǻ, to be compared to the experimental value of 5.44 Ǻ at  K [55]. After relaxation, the lattice constant of GaP is recovered almost everywhere along the direction in the GaP layer, whereas the GaAs layer is strongly strained along the direction with an average lattice parameter  Ǻ. For GaAs bulk,  Ǻ in our DFT simulations, in very good agreement with the experimental value of 5.65 Ǻ, at low temperature [55]. It should be noted that fairly corresponds to the value estimated from the linear elasticity theory using experimental elastic constants ( GPa and  GPa).

DFT offers the possibility to study the structural properties beyond classical elasticity, especially close to the interfaces between materials [11, 15]. Figure 7 represents the variation of the atomic interlayer spacing close to the P-Ga-As interface of the relaxed stack shown in Figure 6. The same stack is simulated using the semiempirical VFF approach [16]. VFF is optimized to describe semiconductor heterostructures at room temperature. Thus, the VFF atomic interlayer spacing values are then shifted to slightly higher values when compared to the DFT results. Both methods show nevertheless that bond stretching undergoes a small fluctuation close to the interface. This fluctuation only extends over roughly one monolayer in the two layers apart from the interface. The layer spacing fluctuation can be interpreted by using an atomic charge analysis. Indeed, the Hirshfeld population analysis [20, 42] is based on the separation of the electron density in real space close to atoms. Conceptually simple, it should be however used with some caution [43]. It was used to describe qualitatively the variation of atomic polarizabilities and atomic charges in III-V semiconductor nanostructures with heteropolar bonds [44]. Figure 8 shows the variation of the atomic charges in short-period GaAs(4)/GaP(4) SLs. Atomic charges calculated for P (−0.15 e) and As (−0.18 e) atoms correspond to the ones obtained in GaP and GaAs bulk, respectively. This also applies for Ga atoms ( e or +0.18 e), except for those located exactly at the interface (+0.16 e). The stretching of the P-Ga and Ga-As bonds connected to the Ga atom at the interface can thus be related to their reduced polarizabilities. This effect is well reproduced by the VFF computations. It should be stressed that advanced empirical force field models are now developed beyond the VFF method in order to get a fine agreement with DFT results for elastic constants and LO-TO phonon splitting [56].

The drop of the interface’s potential between materials in heterostructures is a very important parameter for PV structure’s design. Modern theory of band-lineup in semiconductors relies on the model solid approach of van de Walle and Martin [5759]. The calculation of the valence band-lineup between two semiconductors A (GaAs) and B (GaP) is based on the construction of a superlattice with A/B pseudomorphic interfaces as presented in Figures 6 and 9. The drop at the A/B interface and the planar averaged ab initio potential is then calculated [57]. To complete the computation of valence band discontinuity , one has also to perform the band calculations for bulk and shift the averaged ab initio potentials and valence band maxima, and . The valence band discontinuity is then calculated by transitivity . The conduction band energy levels are represented schematically in Figure 9, but they are not used in the valence band-lineup computation. The valence band-lineup must also be corrected to account for strain and SOC effects [59] or for the valence band self-energy corrections at the G0W0 level [11, 15]. The one-shot G0W0 self-energy correction is known to imperfectly correct the band-gap calculated at the LDA level. The quality of the G0W0 approximation is indeed closely tied to the quality of the LDA eigenfunctions at the starting point. The quasiparticle self-consistent GW (QSGW) method [12] results in state-of-the-art predictions of excited-state properties. However, the additional corrections associated with QSGW self-consistency mostly affect conduction band levels in conventional semiconductors [12].

We apply in the following the model solid approach for the GaAs/GaP interface and would like to point out that a similar LDA + G0W0 study was already presented for the InAs/GaP interface taking into account strain and SOC [15]. The LDA + G0W0 study of unstrained bulk GaAs and GaP (Figure 10) yields as expected a direct band-gap for GaAs and an indirect band-gap for GaP. A detailed inspection of the G0W0 corrections (Table 1) demonstrates that they are required to recover the right order for the X and L band-gaps in GaP. In addition, the effects are more pronounced for the two compounds at Γ than at X. This point is very important for GaAsP(N)/GaP(N) SL [53] and (In)GaAs/GaP quantum dots [60], where the observed optical transitions result from a fine competition between the Γ and X conduction band states. The DFT simulations carried out only at the LDA level for superlattices are then not expected to yield accurate predictions especially for the conduction subbands. As DFT + LDA fails to reproduce correctly the electron effective masses [12], the quantum confinement is not indeed reproduced correctly for small nanostructures.

LDA + G0W0 studies have been recently applied directly to very short-period SL instead of using the corrections for bulk materials [14, 61]. Major error in that case however lies in the assumption that the self-energy in the central part of the superlattice can already be identified with the self-energy of the bulk. The GW self-energy contains indeed the huge image effect due to the dielectric discontinuities at the interfaces of short-period SL. This can be partially corrected using a semiempirical approach based on the dielectric constant discontinuity [61]. In this work, we present LDA + G0W0 studies for bulk GaAs and GaP, but taking into account strain effects. The computed G0W0 valence band shifts are equal to −0.08 eV and −0.25 eV for unstrained GaAs and GaP, respectively. The interplay of strain and many-body effects can be investigated for bulk materials. LDA + G0W0 calculations were performed for a biaxially strained GaAs on GaP (001) (Figure 10). In order to understand the effect of the strain state, the strain tensor was decomposed into a hydrostatic component and a pure biaxial component . The energy shifts and definitions of deformation potentials are given in Figure 11. In addition, the interpretation of the electronic band diagram shown in Figure 10 must take into account the band folding effect associated with the doubling of the primitive cell from the rhombohedral system to the tetragonal one. Electronic energy band-gaps deduced for strained GaAs from LDAG0W0 at the and Γ points are almost equal. Calculations based on experimental deformation potentials predict an indirect band-gap at point (Table 2, [53, 62]). The G0W0 valence band shift for GaAs is equal to −0.19 eV.

The planar averaged and fully averaged ab initio potentials are calculated independently in each crystal bulk, for example, and , and for the GaAs(8)/GaP(8) SL (Figure 9, [20, 64]). The shifted planar averaged potentials for bulk GaP and strained bulk GaAs are shown in Figure 12, respectively, in green and red. The potential curve of GaP has been downshifted by −0.66 eV to fit into the respective part far from interface of the potential curve of the GaAs(8)/GaP(8) SL. The potential curve of strained GaAs has been upshifted by 0.62 eV. The effective lineup of the LDA averaged potentials for GaP and strained GaAs can then be deduced: . Taking into account the energy shifts for bulk and , the valence band-lineup in a GaAs/GaP heterostructure is finally computed at the LDA level:

Adding the valence band G0W0 corrections, one finally ends up with a valence band-lineup value of . The GaAs/GaP valence band-lineup calculated from experiment [55] including strain effects and without SOC is equal to 0.55 eV. Another study using a similar method at the LDA level with a somewhat smaller plane wave cut-off (220 eV) yields a comparative value of 0.61 eV [65]. As three independent DFT calculations are needed in the model solid approach (Figure 9), such an agreement is satisfying.

Finally, it is possible to calculate at the LDA level the valence band spin-orbit splitting of GaAs and GaP using pseudopotentials constructed from the HGH scheme [25]. Theoretical results are in good agreement with experiment [55] for GaAs 0.35 eV (exp. 0.34 eV) and GaP 0.04 eV (exp. 0.08 eV). State-of-the-art DFT calculations have been recently developed in order to fully take into account SOC in the single-particle Green function G and the screened interaction W of the GW approximation [66]. It is possible in such a way to correctly predict the reverse band ordering of II-VI semiconductor alloys containing Hg atoms. These compounds are important for the design of topological insulators [67, 68].

4. Conclusions

We have demonstrated that DFT simulations yield relevant information on PV materials. Quantitative design of PV devices is beyond the possibility of nowadays DFT codes but greatly benefits from reliable physical predictions of structural and optoelectronic properties of bulk materials. Studies of small heterostructures may give also additional information on atomistic relaxations, interfaces, and electronic band-lineups. DFT methods are efficiently used in connection with atomistic simulation tools, like the VFF and tight binding codes. Much effort has been done recently by the DFT community to take into account precisely many-body effects for bulk and small nanostructures or to develop a density-functional perturbation theory. These improvements are also promising for the future of PV materials.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgments

This research was performed using HPC resources from GENCI CINES and IDRIS 2013 (c2013096724). It is also supported by Agence Nationale pour la Recherche (PEROCAI Project ANR-2010-04 and MENHIRS Project ANR (Grant no. ANR-2011-PRGE-007-01)).