Pingping Jiang

Marie-Christine Record

Pascal Boulet

First-principle calculations on CuInSe 2 /AlP heterostructures

Keywords: Buffer layer, heterostructure, band alignment, absolute deformation potential, biaxial strain, topological properties

Heterostructures based on CuInSe 2 absorber with AlP buffer have been modeled for the first time by different stacking schemes and interfacial terminations. Mechanical, electronic and topological properties of CuInSe 2 /AlP heterostructures along [001] direction were investigated using full potential linear augmented plane wave (FP-LAPW) method. Optimal interfacial distance, thermodynamic stability and band alignment have been computed systematically. Heterointerface at "on-top" stacking is the most thermodynamically stable compared with other stackings. The magnitudes of core level binding energy differences are in accordance with interfacial adhesion energies. Absolute deformation potentials ADPs of band energies as a function of (001) biaxial strain were determined. With redressing ADPs, the revised valence and conduction band offsets have increased and decreased, respectively. AlP has strong resistance to biaxial strain in the evolution of band energy edges. Topological properties and density of states were used to investigate the relationship between interatomic interactions and energy band evolution. The most stable CuInSe 2 /AlP heterostructures have "spike-like" band offsets, with the one terminated as Se-Al characterized by a flat conduction band offset of 12 meV. This comes from the strong "Cu d-Se p" states hybridization. The theoretical maximum conversion efficiency of AlP-based "absorberbuffer" heterostructure is 27. 39% at 0.5 μm thickness.

Introduction

Thin-film heterojunction solar cells have become alternatives to traditional silicon-based ones due to superior properties. 1 CuInSe 2 is seen as a promising solar cell candidate with direct bandgap, high absorption and reduced component toxicity compared to other materials (CdTe, PbS, etc). 2 Buffer layer plays an important role in lattice adaptation, adhesion improvement and stress relaxation between itself and the subsequent thin film 3 . A world-record efficiency of 22.3% has been achieved based on Cu(In,Ga)Se 2 absorber and CdS buffer 4 . Alternative buffer materials have been investigated for decades, primarily focusing on photo-generation, Cd replacement with environmentally friendly elements, and thermal stability adjustment. 5 Bandgaps and band alignments are critical parameters for the design of solar cells. Bandgap determines a certain transparency at near-infrared spectrum, hence the degree of light absorption.

Band alignments characterize the energy discontinuities in band edges, which serve as the basis for controlling transport properties and quantum confinement. [START_REF] Sze | Physics of semiconductor devices[END_REF] Buffer materials with satisfied band alignments and well-matched lattice parameters with absorber layers are essential. [START_REF] Hamri | [END_REF] Besides, the effect of epitaxial strain on electronic structures, absolute deformation potentials (ADPs) 8 and band alignments 9 should be focused on as well. Some progresses have been achieved in the calculation of these parameters from the first-principle theory. Vurgaftman et al 10 found that the total hydrostatic deformation potential is proportional to the pressure coefficient of direct bandgap. S.

Kadantsev et al 8 found that the nonzero ADPs of reference core levels have reduced band offsets under compressive biaxial strain.

Figure 1a

presents lattice constants and bandgaps of III-V and II-VI semiconductor compounds, as taken from reference. [START_REF] Sze | Physics of semiconductor devices[END_REF] Considering bandgaps and lattice matches, cost and environmental friendliness, "Cd-free" buffer material, AlP, has been suggested as suitable buffer material for its large bandgap and small lattice misfit with CuInSe 2 (less than 6%). Besides, AlP has higher refractive index than CdS, i.e. greater transmittance, in the range of [1.0;4.5] eV (Figure 1b), which contributes to the absorption degree of absorber layer. Meanwhile, when the photon energy hf is below 4.5 eV, CuInSe2 has a higher absorption coefficient; when the photon energy hf is above 4.5 eV, AlP has a higher absorption coefficient (Figure 1c). The well-combined spectrum absorption of CuInSe 2 /AlP heterostructure maximizes its whole optical ability. Until now, by contrast to other III-V compounds, AlP has not been studied as a buffer layer. In order to provide a foundation for its application, a comprehensive investigation of the properties of CuInSe 2 /AlP heterostructure needs to be done.

In this work, heterostructures based on CuInSe 2 absorber with AlP buffer layer were modeled along [001] direction by using different stacking schemes and interfacial terminations. This direction was chosen as it corresponds to the smallest lattice mismatch (5.6%) between the two layers.

The interfacial distance d 0 and work of adhesion W ad have been quantitatively examined to acquire a thermodynamically stable interface. Band alignments and ADPs of band energies as a function of (001) biaxial strain have been systematically computed. Recently, topological analyses of electron density based on the quantum theory of atoms in molecule (QTAIM) [START_REF] Bader | Atoms in molecules: a quantum theory[END_REF] have been used for interpreting interatomic interactions in real space for photovoltaic materials. [START_REF] Jiang | [END_REF] In combination with density of states, the stability and growth mechanism of heterostructure have been analyzed.

Ultimately, the optical responses of heterostructures based on AlP and CdS buffer layers, as well as that of pure absorber have been estimated and compared as film thickness varies.

Computational details

First-principle calculations were carried out by full potential linear augmented plane wave method (FP-LAPW) as implemented in program WIEN2k. [START_REF] Blaha | WIEN2k, an augmented plane wave+local orbitals program for calculating crystal properties Karlheinz Schwarz[END_REF] The WC-GGA 14 functional has been used to describe exchange-correlation energy in structure optimizations. The convergence criteria on total energy, charge and force have been set to 10 -5 Ry, 10 -4 and 2 mRy/bohr, respectively. In order to improve the bandgap, which is underestimated with GGA functional, the modified Becke-Johnson exchange potential (mBJ) has been used. 15 The planewave cutoff parameter R MT K max in the interstitial region was set to 7. A mesh of 1500 kpoints determined using the Monhkorst-Pack scheme was used to sample the Brillouin zone. modeled using slabs composed of (CuInSe 2 ) m /(AlP) n with a thickness of 15.0 Å vacuum layer. This vacuum layer was applied to separate free surfaces, avoiding interactions between periodic images.

AlP and CuInSe

The lattice constant m and n correspond to the numbers of atomic layers in CuInSe 2 and AlP slabs, respectively. To meet the boundary condition and lattice coherence, the in-plane lattice constants (a and b) of heterostructure were both set to the mean value of bulk CuInSe 2 and AlP lattice constants. In such way, AlP side and CuInSe 2 side are under tensile and compressive strain respectively. The dependences of band edges on biaxial strain have been simulated as deformation ratios vary from -4% to 4% with respect to the equilibrium in-plane lattice constants a and b (with a=b) of CuInSe 2 and AlP, respectively. The chemical potential of each species in our model can be obtained by the first derivative of internal energy to particle number by neglecting the entropy and volume changes. Topological properties of electron density calculated in WIEN2k have been studied by program CRITIC2 16 . Spectroscopic limited maximum efficiency (SLME) has been valued as a useful indicator to characterize photovoltaic abilities. 17 In solar cell applications, conversion efficiency depends on the thickness of absorber layer. 18 The SLME of AlP-and CdSbased heterostructures, as well as pure CuInSe2 thin film, have been calculated under AM 1.5G solar flux.

Surface investigations

The optimized lattice parameters (a and c), bulk moduli (B), bandgaps (Eg) and elastic constants (C ij ) of CuInSe 2 and AlP, are listed in Table 1, together with the corresponding data taken from literature. The calculated lattice constants, bulk moduli and axial elastic constants with WC-GGA functional are well-matched to the literature data, whereas the shear modulus (C 44 and C 66 )

of CuInSe 2 and the bandgaps of both CuInSe 2 and AlP are overestimated and underestimated, respectively. It is found that the bandgaps of CuInSe 2 and AlP have been improved by using the mBJ potential, which agree well with the experimental data. This validates our computational settings and the mBJ potential will be used hereafter. AlP(001) slab with symmetric atomic layers were built to eliminate the spurious dipole effects.

After full relaxation, the relative deviations of interlayer spacing with respect to bulk (Δd ij ) were calculated as N varies, as listed in Table S1 (refer to the supplement information Table S1). The

Δd ij values at the middle layers have been used to determine the slab convergence, as plotted in the Δd ij values tend to reach a plateau around 1% for two terminations. This 13-layers AlP(001) slab is thick enough to exhibit bulk-like interior, its surface energy σ can be written as:

𝜎 = ! "# ! #𝐸 $%&' -∑ 𝑁 ( 𝜇 ( $%&' ()#%,+ -𝑇𝑆 + 𝑃𝑉. (1)
where, 𝐸 $%&' is the slab energy, N x and 𝜇 ( $%&' are the numbers and chemical potentials of atoms in slab, respectively, and A s is the surface area. P, V, T and S are the pressure, volume, temperature and entropy, respectively. Under atmospheric pressure, the products of TS and PV are zero at 0 K.

As suggested in Ref. 25 the slab surface is in equilibrium with bulk, which means 𝜇 #%+ ',%-= 𝜇 #% $%&' + 𝜇 + $%&' . Thus, equation (1) can be rewritten as:

𝜎 = ! "# ! #𝐸 $%&' + (𝑁 + -𝑁 #% )𝜇 #% $%&' -𝑁 + 𝜇 #%+ ',%-. (2)
With the aim to prevent the slab from breaking into energetically more stable substances, the chemical potential of atoms in the stabilized slab 𝜇 $%&' should be lower than those in bulk phases 𝜇 ',%-. Following thermodynamics rules, 𝜇 #%+ ',%-equals to the sum of chemical potentials of atoms in bulk and an energetic term corresponding to compound formation energy ( ∆𝐺 . / ), reading

𝜇 #%+ ',%-= 𝜇 #% ',%-+ 𝜇 + ',%-+ ∆𝐺 . / .
The formation energy ∆𝐺 . / (𝐴𝑙𝑃) is approximated as the energy difference per formula ∆𝐸 . / (𝐴𝑙𝑃) between Al, P and AlP. Thus, the chemical potential difference Similarly, the convergence tests of CuInSe 2 (001) slabs, including Cu/In termination and Se termination (inset of Figure 3a), were conducted as well. After full relaxation, the interlayer spacing relative deviations with respect to bulk cannot be extracted precisely due to the great internal atomic displacements. Aiming to rapidly and accurately obtain the appropriate number of atomic layers, the incremental energy approach 27 was applied as the number of atomic layers N increases. The energy difference 𝐸#𝑁 0 . -𝐸(𝑁 1 ) versus the increment of atomic layers 𝑁 0 -𝑁 1 was calculated, as plotted in Figure 3a. The energy difference first decreases and then keeps steady within an error range of 0.02 Ry. When atomic layers number reaches 13, the corresponding energy differences between 13-layers and 9-layers of Cu/In-and Se-terminated slabs are the lowest as well as the closest to the calculated bulk energy. Hence, a 13-layers slab has been considered thick enough to maintain a bulk-like interior. Note that the numbers of Cu atoms and In atoms are equal in the slab, i.e. 𝑁 2, = 𝑁 34 . Thus, its surface energy σ can be written as:

(𝜇 #% $%&' -𝜇 #% ',%-) is
𝜎 = ! "# ! :𝐸 $%&' + : ! " 𝑁 56 -𝑁 2, ; × #𝜇 2, $%&' + 𝜇 34 $%&' . - ! " 𝑁 56 𝜇 2,3456 " ',%- ; (3) b a
Similarly, the 0 K formation energy ∆𝐸 .

/ of CuInSe 2 can be obtained through the energy difference per formula between Cu, In, Se and CuInSe 2 , which is -1.80 eV for two terminations.

The chemical difference of (𝜇 2, $%&' -𝜇 2, ',%-+ 𝜇 34 $%&' -𝜇 34 ',%-) is in the range shows CuInSe 2 unit cell.

Interface investigations

Interface modeling and surface adhesion

The CuInSe 2 /AlP heterostructures with two types of interfacial bondings, i.e. Se-Al and Cu/In-P, were modeled based on the above convergence tests. In effect, since the adjacent layers share neither common anion nor cation, these interfacial bondings unnecessarily lead to the same band offsets. 28 Three stacking schemes for each interfacial termination, designated as "stack1", "stack2" and "stack3", are under consideration. As presented in Figure 4, "stack1", "stack2" and "stack3" means that atoms from the first atomic layer of CuInSe 2 (001) slab sit at the "diagonal" midpoint, b a at the "lateral" midpoint and directly "on-top" of atoms from the first atomic layer of AlP(001) slab, respectively. Hereafter, the (001) plane index of CuInSe 2 and AlP will be omitted for simplification. The work of adhesion, W ad , as characterized by the energy difference between separate slabs and their combined heterostructure, has been used to portray the heterointerface properties, such as strength and stability. The W ad of CuInSe 2 /AlP is given by, 25 𝑊 &7 (𝐶𝑢𝐼𝑛𝑆𝑒 " 𝐴𝑙𝑃 ⁄ ) = #∑ 𝐸 is the heterostructure energy and A i is the interfacial area. and Cu/In-P-stack1 is observed in comparison with that of -stack2 and -stack3, for which the optimal d 0 and corresponding maximum W ad are both well-matched between the two terminations.

At the optimal interfacial distance, heterostructures have been relaxed until the total residual forces on each atom were less than 2 mRy/Bohr. For the sake of time-and cost-saving, only the b d c a four atomic layers near interface (enclosed within the red dash lines in Figure 4a and4b) have been considered. As shown in Table 2, the work of adhesion after relaxation still follows the same order as the one generated from UBER. As discussed before, the maximum surface energies of Aland P-terminated AlP slabs can reach up to 1.37 J/m 2 and 1.50 J/m 2 , respectively. Correspondingly, their work of adhesion W ad (AlP/AlP) along [001] direction can be obtained according to equation (4), i.e. 2.74 J/m 2 and 3.00 J/m 2 respectively, which are twice the surface energies since the interfacial energy of the "AlP/AlP heterostructure" is viewed as zero. In the thermodynamic system, substrate with sufficient interfacial stability is indispensable for heterogeneous nucleation.

Compared to W ad (AlP/AlP), the W ad of Se-Al-stack3 and Cu/In-P-stack3 heterostructures are larger, which means that CuInSe 2 thin film is more likely to form on AlP substrate than AlP itself.

Band alignment

The valence band offset ∆𝐸 > is calculated following the same procedure as that in core-level photoemission measurement, 30 which is defined as:

∆𝐸 > (𝐶𝑢𝐼𝑛𝑆𝑒 " 𝐴𝑙𝑃 ⁄ ) = ∆𝐸 ?@A,2 2,3456 " -∆𝐸 ?@A,2 # #%+ + ∆𝐸 2,2 # (5)
where, ∆𝐸 ?@A,2 During the electron transition process in photovoltaic devices, effective transportation should happen before the excited electrons recombine or fall back into a lower energy level. CuInSe 2 /AlP interface is characterized by high defect formation mainly due to the lattice and thermal expansion mismatches 31,32 . These defects play a fundamental role in solar cell performance. The heterointerface should exhibit minimum recombination and current loss thus allowing for the maximum open-circuit voltage. According to Courel et al. 33 , interfacial recombination depends on the combining velocity and minority carriers concentration, defined as follow:

𝐽 14;6< = 𝑞𝑆𝑛 1D2,3456 " #1 + 𝛾 EF5 𝑒 D(∆I ( J∆I ) ) "-L ⁄ . (6)
where, 𝐽 14;6< is the reverse saturation current density due to recombination, q is electron charge, S is recombination velocity, 𝛾 EF5 is density of states enhancement factor between CuInSe 2 and AlP.

So, the influence of band offsets on minority carrier recombination can be demonstrated.

Indeed, "spike-like" heterojunction (ΔE c >0) has less interfacial recombination than "cliff-like" one (ΔE c <0). Interface owning positive CBO bears more favorable carrier transportation properties than that owning negative one. However, large CBO is viewed as a sufficiently high barrier so that majority of electrons generated in the absorber layer will recombine before passing through. This electron-hole recombination at the absorber surface is detrimental to open-circuit voltage. Besides, small VBO causes significant electron confinement without hindering the hole injection in a lightemitting device. 34 This will increase the chance of recombination as well. Thus, band offsets should be strictly limited in a certain range. Though Cu/In-P-stack3 has slightly larger interfacial stability than Se-Al-stack3, its large CBO and small VBO greatly deteriorate the performances regarding band lineup and external quantum efficiency.

Absolute deformation potentials (ADPs) under biaxial strain

In the above calculations, there is a basic assumption that the absolute core level deformation potential is zero, therefore the band alignment between solids is intrinsic. However, there is a great uncertainty on deformation potentials, which comes from the strain effects in pseudo-morphically grown layers. 8 To precisely predict band offsets, the ADPs of energy bands, accounting for band edges and core levels, have been considered in the following calculations. The average in-plane lattice constant at interface results in compressive strain at CuInSe 2 slab and tensile strain at AlP one. For a given deformation ratio, their optimized out-of-plane lattice constants and volumes were obtained (refer to supplement information Figure S1). As proposed in Ref. 8 , the change in band energy under biaxial strain is defined as 𝑎 ? '&47 . The ADP of VBM can be expressed as 𝑎 ? ?@A = 𝑎 ? ?@A/!$ + 𝑎 ? !$ , where 𝑎 ? ?@A/!$ is the deformation potential of VBM with respect to E ref , and 𝑎 ? !$ is the reference core level and was determined in the previous natural band offsets calculation.

Figure 6 shows the biaxial strain-modified band energies, including CBM, VBM and core levels, of bulk AlP and CuInSe 2 . As strain goes from compressive to tensile, the band energies decrease continuously, which gives negative ADPs. The nonnegligible ADPs, when subjected to biaxial strains, emphasize the importance of controlling lattice mismatch between thin film and substrate. By fitting band energies to the polynomial of third-degree, the ADPs can be determined from the first derivative of the fitted curves with respect to volume. The 𝑎 ? ?@A , 𝑎 ? 2@A of AlP and CuInSe 2 at the average in-plane lattice constant (marked as "*" in Figure 6) are -82.4 meV/Å 3 , -100.9 meV/Å 3 and -39.2 meV/Å 3 , -20.4 meV/Å 3 , respectively. Taking the core level deformation potential into account, the 𝑎 ? ?@A/!$ of AlP and CuInSe 2 are -0.3 meV/Å 3 and -10.2 meV/Å 3 , respectively. AlP has an equivalent variation of core level and valence band edge when subjected to strain. This strong resistance against strain contributes to its potentiality in buffer layer applications. By contrast, CuInSe 2 has an unbalanced variation of core level and valence band edge. At the average lattice constant, the absolute band energies can be given by the product of ADPs and volume deformation. Associating the biaxial strain modified CBM and VBM with the reference core levels and equation ( 5), the revised VBO ∆𝐸 > <6 and CBO ∆𝐸 B <6 can be expressed as:

a b ∆𝐸 /(0) 1+ = ∆𝐸 /(0) + 𝑎 " "(&)#$ d𝑉(CuInSe ! ) -𝑎 " "(&)#$ d𝑉(AlP) (7)
where, ∆𝐸 >(B) is the natural valence and conduction band offsets obtained previously without considering band energy deformation. The revised band offsets are listed in Table 3, marked by "re". Regarding Se-Al termination, for a given raw band offset, the revised VBO has increased by 94.0 meV, as a result of the energy upshift of VBM towards core level in CuInSe 2 which amounts to 93.1 meV, and the negligible energy downshift of -0.9 meV in AlP. The revised CBO has decreased by 102.3 meV due to the large downshift of CBM towards VBM in CuInSe 2 , which is -171.4 meV, coupled with a small downshift of -69.1 meV in AlP. The revised VBO and CBO for Cu/In-P heterostructure termination can be obtained as well as soon as the energy shifts between VBM, CBM and core level are determined. as suggested in Ref. 35,36 . Accordingly, the considerable "spike-like" CBO of Cu/In-P-stack3 shall constrain its interest in heterostructure assessment, since CBO is recognized as the depth of quasiparticle confinement. As mentioned in Ref. 37 , the discrepancy in strain states between two neighboring heterointerfaces is directly correlated to the asymmetric charge distribution feature.

Electrons released from the CuInSe 2 absorber layer are supposed to cross the barrier in the short term before getting recombined or degenerated. In view of that, Se-Al-stack3 heterostructure has become our best choice due to its favorable band lineup. Furthermore, its thorough investigations on strain modified band edges provide us perspectives on choosing and determining appropriate buffer layer in thin film manufacturing. These strain-mediated evolutions of band edges could offer additional degrees of freedom for achieving the desired functionality and utility.

Spectroscopic limited maximum efficiency SLME

According to Crovetto and Hansen 42 , low-efficiency cells have large "cliff-like" offsets, whereas most high-efficiency cells have "spike-like" or nearly flat offsets. According to the preceding discussions, the "on-top" stacking terminated as Se-Al at interface is foreseen as the most promising heterostructure. To better characterize both the designated "absorber-buffer" heterojunction and the pure CuInSe2 absorber their photovoltaic responses have been calculated.

In addition to the AlP-based heterostructure, a CdS-based one has been modeled with the same interfacial condition for the purpose of comparison. By computing the absorption coefficient in and CdS serve as buffer layers. As seen from Figure 9b, the maximum efficiencies of heterostructures and pure CuInSe2 increase rapidly at the early stage and then tend towards the converged value. It is found that 0.3 μm CuInSe2 absorber is seen as thick enough to attain a converged efficiency of 27.42%, which is close to the results of Ref. 43,44 under the same condition.

Similarly, 0.5 μm Se-Al-stack3 and Cd-Al-stack3 heterostructures have attained their converged efficiencies, 27.39% and 27.41%, respectively. These two heterostructures possess not only high (close to pure CuInSe2 absorber) but also similar optical properties. In other words, AlP can be utilized as a promising alternative in buffer material in pursuing high efficiency and low-cost solar cell. Therefore, the addition of AlP buffer merely affects the whole scale of optical properties, which undoubtfully, benefits from the well-matched interface condition, as illustrated before.

Author information

Corresponding Author

*E-mail: m-c.record@univ-amu.fr.

Acknowledgments

This work is financially supported by the China Scholarship Council (CSC). This work was granted access to the HPC resources of the Centre Informatique National de l'Enseignement Supérieur (CINES), Montpellier, France under allocation A0050806881 made by the Grand Equipement National de Calcul Intensif (GENCI). It was also granted access to the HPC resources of Aix-Marseille Université financed by the project Equip@Meso (ANR-10-EQPX-29-01) of the program "Investissements d'Avenir" supervised by the Agence Nationale de la Recherche.

Figure 1 .

 1 Figure 1. (a) Bandgaps and lattice constants of II-VI (red color), III-V (blue color) and CuInSe 2 (magenta) at ~2K, as taken from Ref. 6 . Abbreviations of "zb", "wz" and "rs" indicate zincblende, wurtzite, and rock salt crystal structures, respectively. Compounds without abbreviation crystallize in zb structure or related one (CIS). (b) Calculated refractive index of AlP and CdS and (c) absorption coefficient of AlP and CuInSe 2 .

  2 crystallize in cubic (space group F-43m) and tetragonal structures (space group I-42d), respectively, as shown in the inset of Figure 2a and 3a. The heterostructures wereb a c

Figure 2a .

 2a Figure 2a. For low numbers of layers, Δd ij decreases sharply as N increases. When N attains 13,

  in the range [∆𝐸 . / ; 0], being [-1.18;0] eV and [-1.21;0] eV for Al-and Ptermination respectively.The σ evolution of AlP(001) slab as a function of chemical difference (𝜇 #% $%&' -𝜇 #% ',%-) can be obtained from equation (2), as soon as the chemical potentials of Al and AlP in bulk phases are determined. As plotted in Figure2b, the σ values of Al-and P-termination are in the ranges 1.06-1.37 J/m 2 and 1.17-1.50 J/m 2 , respectively. Besides, there is a cross-point between them where chemical difference is -0.82 eV. Lower surface energy means larger stability.26 In the ranges of [-1.21;-0.82] eV and [-0.82;0] eV, P-and Al-termination show opposite behaviors in evaluating the preferential surface stability.

Figure 2 .

 2 Figure 2. Relative deviations of interlayer spacing with respect to the number of atomic layers in AlP(001) slab for Al-and P-terminations (a) and their corresponding surface energies at 13-layers as a function of chemical difference (𝜇 #% $%&' -𝜇 #% ',%-) (b). The inset of (a) shows AlP unit cell.

  [2∆𝐸 . / ; 0] . The σ evolution as a function of chemical difference can be deduced from equation (3), as soon as the chemical potentials of Cu, In and CuInSe 2 in bulk phases are determined. As plotted in Figure 3b, the σ values of Cu/In-and Se-termination are in the ranges 0.84-1.70 J/m 2 and 0.16-1.02 J/m 2 , respectively. There is also a cross-point where chemical difference is -0.38 eV. Se-termination and Cu/In-termination behave reversely in seeking for favorable stability in the ranges of [-3.60;-0.38] eV and [-0.38;0] eV. Therefore, both Al-and P-terminated AlP surfaces together with Cu/In-and Se-terminated CuInSe 2 surfaces need to be considered in the forthcoming heterostructure modeling.

Figure 3 .

 3 Figure 3. Energy difference 𝐸#𝑁 0 . -𝐸(𝑁 1 ) versus the increment of atomic layers number 𝑁 0 -𝑁 1 in CuInSe 2 (001) slab for Cu/In-and P-terminations (a) and their corresponding surface energies at 13-layers as a function of chemical difference (𝜇 2,$%&' -𝜇 2, ',%-+ 𝜇 34 $%&' -𝜇 34 ',%-) (b). The inset of (a)

Figure 4 .

 4 Figure 4. Stacking schemes of CuInSe 2 /AlP heterostructures: (a) Se-Al termination; (b) Cu/In-P termination. "stack1", "stack2", and "stack3" represent "diagonal" midpoint, "lateral" midpoint, and "on-top" stacking as depicted in (c) for Se-Al termination and (d) for Cu/In-P one.

Figure 5 .

 5 Figure 5. Work of adhesion W ad vs. interfacial distance d 0 of Se-Al (solid symbols) and Cu/In-P (hollow symbols) terminations for each stacking. Solid and dash lines are their smoothing curves.Color code : stack 1 in red, stack 2 in black, stack 3 in purple.

Figure 6 .

 6 Figure 6. Evolution of band energies, for VBM, CBM and core levels of bulk AlP (a) and CuInSe 2 (b) as biaxial strain ratio varies. The star symbol "*" shows the average lattice constant ratio.

Figure 7

 7 Figure 7 depicts the schematic diagram of revised band offsets. Se-Al-stack1 heterostructure has been changed from type-I into type-II due to ADP correction. Se-Al-stack3 and Cu/In-P-stack3, as the most stable heterostructures, have revised CBOs of 12 meV and 817 meV respectively. The nearly flat "spike-like" CBO of Se-Al-stack3 will boost open-circuit voltage and electron transition,

Figure 7 .

 7 Figure 7. Schematic diagram of revised band offsets and corresponding types of Se-Al termination (a) and Cu/In-P termination (b) of CuInSe 2 /AlP heterostructures, and illustration of the migration of electrons and holes at the type-I and type-II heterointerfaces (c).

Figure 8 .

 8 Figure 8. Schematic diagrams of Se-Al-stack3 (a) and Cu/In-P-stack3 (e) near the interface, and their electron density ρ (b and f), electron density Laplacian ∇ " 𝜌 (c and g), and bond degree (BD=H BCP /ρ BCP ) vs. |V BCP |/G BCP (d and h) at BCPs in the plane (-110). The asterisk sign (*) stands for the longer bond length for pairwise Al-P. The insets in (d) and (h) show total and partial density of states of atoms near interface.

[Figure 9

 9 Figure 9.

Figure 9 .Figure S2 .

 9S2 Figure 9. Calculated Jsc-Voc characteristics (a) and spectroscopic limited maximum efficiency (b) with respect to the thickness of Se-Al-stack3 and Cd-Al-stack3 heterostructures and pure CuInSe 2 under AM 1.5G, namely 100 mW/cm 2 illumination. The range of light absorption wavelength is 280 nm-1200 nm.

Table 1 .

 1 Calculated

			lattice parameters (a and c), bandgaps (Eg), bulk moduli (B), and elastic
	constants (C ij ) of CuInSe 2 and AlP, together with corresponding data taken from literature.
	Bulk		a (Å)	c (Å)	Eg(eV)	B (GPa)	Cij (GPa) C11 C12 C13 C33 C44 C66
	CuInSe2	This work	5.786 11.606	0.22 WC 1.06 mBJ	73.7	94	61 67	83 109 118
		Calc.	5.782 a 11.622 a	-	-	71 b 45 b 45 b 63 b 46 b 47 b

b Reference 20 . c Reference 21 . d Reference 22 . e Reference 23 . f Reference 24 . g Reference 10 .

Table 2 .

 2 For both two terminations, the maximum values of W ad for different stacking schemes follow the order W ad (stack3)>W ad (stack1)>W ad (stack2). "On-top" stacking at the optimal interfacial distance is the most thermodynamically stable, with the highest W ad of 4.405 J/m 2 and 4.467 J/m 2 for Se-Al and Cu/In-P terminations respectively. Besides, a considerable difference in W ad between Se-Al-stack1

Table 2 .

 2 Optimal interfacial distance d 0 with corresponding W ad obtained from smoothing curves in Figure 5, as well as the W ad after relaxation.

	Interfacial terminations	Stacking schemes	d0 (Å)	UBER Wad (J/m 2 )	Relaxed_Wad (J/m 2 )
		stack1	2.303	4.131	3.816
	Se-Al	stack2	3.081	3.761	3.756
		stack3	2.190	4.405	4.408
		stack1	2.444	3.854	3.650
	Cu/In-P	stack2	3.003	3.837	3.479
		stack3	2.265	4.467	4.629

  2,3456 " and ∆𝐸 ?@A,2 # #%+ are energy separations between core level and valence band maximum of bulk CuInSe 2 and AlP, respectively, and ∆𝐸 2,2 # is the difference in core level binding energies between cations and anions at interface. The average core level of cations and anions at interface is chosen as the reference core level, E ref . For Se-Al interfacial termination, Al 1s and Se 1s are assigned as cation and anion core states respectively. For Cu/In-P interfacial termination, 1/2(Cu 1s+In 1s) and P 1s are set as cation and anion core states respectively. Besides, the conduction band offset ∆𝐸 B can be obtained by ∆𝐸 B = ∆𝐸 C + ∆𝐸 > , where ∆𝐸 C is the calculated bandgap difference between bulk CuInSe 2 and AlP.

Table 3

 3 lists the energy separations from E ref to VBM, core level binding energy differences ∆𝐸 2,2 # , natural band offsets and corresponding types of CuInSe 2 /AlP heterostructures. Core level binding energy difference is correlated with band alignment; the greater value of ∆𝐸 2,2 ∆𝐸 $ %& , ∆𝐸 ' %& and type re represent the revised valence and conduction band offsets and corresponding types of heterostructures when considering the ADPs of band energies under biaxial strain. These results will be discussed in the next section.

	𝐶𝑢𝐼𝑛𝑆𝑒 ! /𝐴𝑙𝑃 heterostructures	Se-Al termination stack1 stack2 stack3	stack1	Cu/In-P termination stack2	stack3
	∆𝐸 "#$,& &'()*+ !				918.2201			1344.0435
	∆𝐸 "#$,& " ,-.				110.5463			152.5341
	∆𝐸 &,& "			-807.5684 -807.5554 -807.5742 -1191.4516 -1191.4119 -1191.4683
	∆E / (𝐶𝑢𝐼𝑛𝑆𝑒 ! 𝐴𝑙𝑃 ⁄	)	1.434	1.608	1.355	0.786	1.326	0.559
	∆𝐸 0 (𝐶𝑢𝐼𝑛𝑆𝑒 ! 𝐴𝑙𝑃 ⁄	)	0.035	-0.139	0.114	0.683	0.143	0.910
	Type			I	II	I	I	I	I
	∆𝐸 / 1+ (𝐶𝑢𝐼𝑛𝑆𝑒 ! 𝐴𝑙𝑃 ⁄ ∆𝐸 0 1+ (𝐶𝑢𝐼𝑛𝑆𝑒 ! 𝐴𝑙𝑃 ⁄	) )	1.528 -0.067	1.702 -0.241	1.449 0.012	0.920 0.590	1.460 0.050	0.693 0.817
	Type re			II	II	I	I	I	I

# , the larger value of ∆𝐸 > . Se-Al-stack1 and -stack3 heterostructures have type-I ("spike-like") band offsets, where both CBM and VBM belong to CuInSe 2 ; Se-Al-stack2 has type-II ("cliff-like") band offset, where CBM and VBM belong to AlP and CuInSe 2 , respectively. Besides, Cu/In-P-stack1, -stack2 and -stack3 heterostructures all have type-I (" spike-like") band offsets. Regardless of the types of terminations, their absolute values of ∆𝐸 2,2 # at interface follow the order stack3>stack1>stack2, which is surprisingly in accordance with their work of adhesions. Table 3. Energy separations (Ry) from core level to VBM of bulk CuInSe 2 and AlP, core level binding energy differences (Ry), band offsets (eV) and corresponding types of CuInSe 2 /AlP heterostructures. *

Conclusions

In this work, we have chosen the AlP compound as a buffer layer for CuInSe 2 thin film. The CuInSe 2 /AlP heterostructures have been modeled by different stacking schemes and interfacial terminations. The optimal interfacial distance, thermodynamic stability and band alignment have been computed. Heterointerface at "on-top" stacking is the most thermodynamically stable compared with other stackings. There is a correlation between band alignments and core level binding energy difference at the heterointerface. The biaxial strain-induced deformation of band energies and hence their absolute deformation potentials have been determined. The revised valence and conduction band offsets of CuInSe 2 /AlP heterostructure have increased and decreased, respectively. Both Se-Al and Cu/In-P "on-top" stackings has "spike-like" band offsets, and the former has flat CBO, 12 meV. Topological analysis of electron density unveils interatomic interactions at heterointerfaces and provides insightful information for material manufacturing. The strong mixing between Cu d and Se p states plays a decisive role in defining the absolute position of valence band maximum and thus the chemical trends. The SLME studies reveal that CuInSe 2 /AlP heterostructure can attain 27.39% of efficiency, which is similar to that of CdS-based heterostructure. This provides a solid basis for AlP as buffer layers in future applications.

Supporting information

The following files are available free of charge.


Pingping Jiang

Pascal Boulet

Marie-Christine Record

A DFT study of the electronic, optical and topological properties of free and biaxially strained CuIn1-xAlxSe2

Keywords: DFT, CuIn1-xAlxSe2, Optical properties, Electron density topological properties, Biaxial strain

The electronic and optical properties of free and biaxially strained CuIn1-xAlxSe2 have been calculated by using the full potential linear augmented plane wave (FP-LAPW) method. Results show that CuIn1-xAlxSe2 has a direct band gap with increasing value as x increases. CuIn0.75Al0.25Se2 has been recognized as the optimal substituted compound in terms of band gap and conversion efficiency. Geometry optimized bond length and bond angle between nuclei terminal atoms (M, Se) and bond critical point, electron density ρ, local energy density, and Laplacian ∇ ! 𝜌 for nonequivalent pairwise M-Se have been examined. High bond deflection with small bond angle gives high possibility of electron transition, which mainly occurs between In and Se orbitals. Under biaxial strains, the optical properties are improved along the lengthened deformation direction.

Introduction

The conversion of solar energy into electricity, namely photovoltaic (PV) effect, has been regarded as a competitive and attractive substitute to other energy conversion forms due to its tremendous advantages 1,2 . However, the photovoltaic conversion efficiency of semiconductors (SCs) needs materials with appropriate band gap (Eg). Hence, due to the materials incapability to absorb photons of energy lower than Eg and to the low conversion efficiency, caused by thermal losses, when photon energy is far higher than Eg, we have to explore more adaptable PV materials 3,4 . The ternary Cu-based chalcopyrite compound, CuInSe2 (CIS), is an interesting material as solar cell absorber layer due to its low cost, high absorption coefficient, excellent optical and electrical properties. Many approaches have been adopted to improve the energy conversion efficiency 5 . However, its narrow band gap and the scarcity and expensiveness of indium constrain its large-scale development. Substituting indium with the abundant and inexpensive aluminum to form the quaternary SC, CuIn1-xAlxSe2 (CIAS), has been considered as a promising alternative with few changes in physical and chemical properties 6 .

A general trend is that the band gap decreases in moving down in a periodic table group 7 . Theoretically, the optical band gap of CIAS can be varied from 1.04 eV to 2.67 eV as aluminum ratio x goes from 0 to 1 8 . Compared to Ga-doped CuInSe2, CIAS allows for reaching higher energy band gap with less aluminum content, therefore fewer lattice distortion and defects. Several experimental researches have been conducted on CIAS, attempting to find the optimal substitute content based on the structural, optical and electronic properties [8][9][10][11] . Aluminum substitution not only increases the band gap, as found by López-García et al 11 on the way to study CIAS (x=0~0.6) thin film, but also changes the carrier density and mobility. As found by Cheng et al 10 , due to low formation energy the main acceptors Cu vacancies decreased as x increased, contributing to such carrier mobility that materials yield the highest photovoltaic performance when x is 0.34. Besides, density functional theory (DFT) has been applied in solar cells area for many years to calculate the structural, electronic and optical properties 12,13 . The agreement between the calculated data and the experimental ones proved the validity of regarding numerical calculations as an adequate method for PV materials studies [14][15][16] . By calculating the optical properties, Tang et al 17 found that CuIn0.75Al0.25Se2 is a suitable PV material in considering both band gap and optical properties.

The quantum theory of atoms in molecule (QTAIM) proposed by Bader et al [START_REF] Bader | Atoms in molecules: a quantum theory[END_REF] has evolved to be an invaluable tool for the chemical interpretation of quantum mechanical data. This theory is built upon quantum mechanical observable quantities, such as the electron density ρ. By inspecting the electron density within the quantum theory framework and analyzing the interatomic interactions through the topological analysis of ρ, a comprehensive characterization from the strong shared-shell (SS) to weak closed-shell (CS) interactions has been concluded. The classification of the SS-and CS-type of interactions is based on the negative and positive sign of the Laplacian of ρ at their bond critical points (BCPs), respectively. However, it was found that such a simple classification in description of these bonds is not meticulous enough, particularly in complex situation [START_REF] Gervasio | [END_REF] . A more specific quantitative scale is then needed. The local energy density obtained by integrating along each basin, the electron density flatness, the charge transfer, and the bond degree have been proposed to better characterize the interatomic interactions [START_REF] Gervasio | [END_REF][20][21] . Saeed et al 22 found that for Ag-based chalcopyrite materials, AgXY2 (X=Al,Ga,In and Y=S,Se,Te), the positive Laplacian of pairwise Ag-Y at BCP characterizes an ionic bond, while the negative Laplacian of the pairwise X-Y at BCP reveals a covalent character. For Cu-based chalcopyrite materials, a comprehensive study of all the topological and energetic properties in the full range of interactions has never been undertaken. Analysis of these properties based on the electron density could help us to clarify the optical behaviors systematically. This can be obtained from the evolution of chemical pairwise interactions as aluminum substitutes indium, giving a chance to shed light on the inner relationship between the topological and optical properties.

In current PV cells, strains originating during their fabrication from the lattice mismatch between the PV materials and the substrates inevitably influence the optical performances. Researches have been conducted on the structure, electronic and optical properties under uniaxial and biaxial strains [23][24][25] . Some researchers even studied the phase stability, electronic properties and chemical bonds of CuIn(Al,Ga)Se2 under high isotropic pressure by DFT 26 . It was found that materials underwent a phase transition that could even lead to nonsemiconductor materials, as with CuGaSe2 for which covalent bonds are strengthened. Studying the influence of strains on photovoltaic materials properties is hence indispensable.

In this work, by DFT calculations we studied the structural, electronic and optical properties of CuIn1-xAlxSe2 (x=0, 0.25, 0.5, 0.75 and 1), and determined the optimal substituting percentage. At this optimal value, we calculated the band gap and optical properties for the material subjected to biaxial strains. Besides, in the aim to unravel the deep relationship between bond interactions and optical properties, a detailed investigation of topological properties based on the electron density has been conducted as strain is applied.

Computational details

We performed first principle density functional theory calculations using the full potential linear augmented plane wave (FP-LAPW) method based on Kohn-Sham equations, as implemented in WIEN2k package [START_REF] Blaha | WIEN2k, an augmented plane wave+local orbitals program for calculating crystal properties Karlheinz Schwarz[END_REF] . The generalized gradient approximation (GGA) approach has been used to describe the exchange-correlation energy [START_REF] Langreth | [END_REF] . CuInSe2 has tetragonal structure (space group I-42d), each selenium ion being coordinated with 2 indium and 2 copper ions, and each indium and copper ion being coordinated with four selenium ions, as shown in the inset of Fig. 1. The core-valence electrons separation has been defined as follows: Cu [Ar]3d 10 4s 1 , In [Kr]4d 10 5s 2 5p 1 , Al [Ne]3s 2 3p 1 , Se [Ar]3d 10 4s 2 4p 4 . The Muffin-tin radii for Cu, In, Al and Se atoms were 2.35, 2.42, 2.24 and 2.24 a.u. respectively. First, we performed structure optimization with PBE-, PBEsol-and WC-GGA. We obtained the optimized lattice parameters by calculating the total energy with respect to the crystal cell volume and c/a ratio in sequence. The fitting result, as shown in Fig. 1, was obtained by using the Birch-Murnaghan equation of state:

𝐸(𝑉) = 𝐸 " + ) 9𝑉 " 16 -. / 𝐵 14703.6 - 6(𝜂 ! -1) # 𝐵 $ + (𝜂 ! -1) ! (6 -4𝜂 ! ) (1) 𝜂 = ) 𝑉 " 𝑉 -. ! " (2)
where Vo, V present the initial and deformed volume respectively, B is the bulk modulus, BP is the derivative of bulk modulus with respect to pressure. After the lattice constants have been optimized, a full relaxation was operated on the whole structure until the total residual forces on all atoms were less than 2 mRy/bohr. Table 1 displays the optimized lattice parameters and bond lengths in comparison to the experimental values of CuInSe2 and CuAlSe2. We can see that, the calculated lattice parameters match the experimental values within a reasonable range of errors and WC-GGA provides the best agreement with respect to experiments. The plane wave cutoff parameter RMTKmax in the interstitial region was set to 7, where RMT is the smallest atomic sphere radii and Kmax is the maximum value of the reciprocal lattice vector in the plane wave expansion. In the interstitial region the charge density and potential were expanded as a Fourier series with wave vectors up to Gmax=14 a.u -1 . The Brillouin zone were sampled with a mesh of 1500 k-points using Monkhorst-Pack grids 29 . These calculating conditions were also adopted in the indium substituted cases. Fig. 2 presents the lattice parameters and their fitting results as Al content x varies. As x increases, lattice a and c decrease linearly, respectively. Bodnar et al 30 studied the CuAlSe2-CuInSe2 phase diagram, confirms that there is a complete solid solution from CuAlSe2 to CuInSe2. Table 2 lists the formation energy of CuIn1-xAlxSe2 when x is 0.25, 0.5 and 0.75 respectively. In-substituted systems are stable due to their negative formation energy.

The biaxial strains have been simulated according to the deformation ratios of ±4%, ±3%, ±2%, ±1%, and 0% with respect to the equilibrium in-plane lattice a perpendicular to c-axis in conditions of constant volume. Negative and positive signs stand for the compressive and tensile in-plane strain, respectively. For each applied constrain the relaxed value of c and u(Se) were calculated. Then, band structure, dielectric constant and absorption coefficient have been calculated as well. As for the topological properties, program CRITIC2 31 has been used for the real space analysis of quantum chemical interactions. In QTAIM, the partitioning in atomic regions, called basin, is determined uniquely by the zero-flux surfaces of the electron density ρ. The kinetic, potential and total energy density, volume and charge can be obtained by the integration in the basins with the electron density calculated by WIEN2k.

Table 2

Formation energy ∇𝑓 of In-substituted system when x equals to 0.25, 0.5 and 0.75, respectively. 3f, the band gap increases continually as x increases. This trend is in accordance with the experiment results 9,17 , though the values have been underestimated by GGA functional because of the inaccurate calculation of exchange-correlation energy [START_REF] Kodan | [END_REF]35 .

Density of states

We use the equilibrium geometries obtained from full relaxation of the CuIn1-xAlxSe2 systems to calculate their total and partial density of states (TDOS and PDOS), which are shown in Fig. 4. The TDOS can be divided into four regions (Fig. 4a). Since the inter-band transition occurs mostly at the vicinity of the Fermi level and no intra-band contributions were added in our study, our main consideration on PDOS focuses in the [-6;10] eV energy region, as shown in Fig. 4b-f. The PDOS presents the contributions to TDOS of different orbitals states for each atom in the system. Fully occupied d orbitals of indium and selenium atoms have not been considered in our investigation. One can see that the valence band is dominated by copper 3d and selenium 4p with the strongest mixing between these orbitals in the upper side of valence band ranging from -2 eV up to the Fermi level. This is consistent with the results in Ref. 14,15 . From -4 eV to-3. characterized by a mixing of In 5s and In 5p, when present in the system, Cu 3d and Se 4p orbitals. Moreover, the densities of these orbitals are at a higher level when x=0.25 than others.

In the aluminum containing compounds, Al 3s and Al 3p orbitals also contribute to the PDOS. The ionic and covalent character of a bond can be determined by the difference of electronegativities between two bonded atoms 36 . Large differences provide an ionic bond due to charge transfers between atoms. By contrast, small differences feature covalent bond due to charge sharing. According to Duffy 37 there is a relationship between band gap and electronegativity, as expressed by 𝜐 %&' = 30000<𝜒 ()* (anion) -𝜒 ()* (metal ion)F , where υmax is the onset of optical absorption, χopt is the optical electronegativity. The decrease of band gap is accompanied by a decrease of ionicity in bonding, which arises from a decreased electronegativity difference. Verma 38 found a linear relationship between band gap and optical electronegativity of solids; CuInSe2 has lower optical electronegativity, 0.588, than CuAlSe2, 0.736, which corresponds to its smaller band gap than the latter one. As the aluminum content x increases, the band gap keeps increasing since the optical electronegativity is enlarged continuously. Thus, for a given incident photon wavelength fewer electrons can jump from valence band to conduction one.

Optical properties

The optical properties of materials, which correspond to the response of materials to an electromagnetic perturbation, can be accessed by the calculation of the frequency-dependent complex dielectric functions 𝜀(𝜔) = 𝜀 + (𝜔) + 𝑖𝜀 ! (𝜔). The imaginary part ε2(ω) is calculated using the following formula:

ε 2 (ω)= e 2 ћ πm 2 ω 2 ∑ ∫ |M cv (k)| 2 δ{ω cv (k)ω} BZ v,c d 3 k (3)
where Mcv(k) are the momentum matrix elements between occupied (valence) and unoccupied (conduction) states, the sums are performed on all the valence and conduction states and the kpoints of the first Brillouin zone. The real part ε1(ω) is calculated using Kramers-Kronig relation:

𝜀 + (𝜔) = 1 + 2 𝜋 𝑃 W 𝜔 8 𝜀 ! (𝜔 8 ) 𝜔 8! -𝜔 ! + 𝑖𝜂 𝑑 9 " 𝜔 8 ( 4
)
where P is the principal value of the integral and η is an infinitesimal number. As ε1(ω) and ε2(ω) strongly depend on the band structure, it must then be properly described. For this, a scissor operator has been used to correct for the underestimated gap with respect to the experimental value. Both the real and imaginary parts of ε(ω) are plotted along the perpendicular (xx) and parallel (zz) directions to the c-axis for CuInSe2 and CuAlSe2 in Fig. 5. The calculated results are compared to experimental ones obtained in Ref. 39,40 from spectroscopic ellipsometry. As it shows, our calculation fits the experimental trends quite nicely, better than Ref. 17 (not shown here), supporting again that our calculation settings lead to an accurate description of the optical properties. Fig. 6 shows the dielectric functions ε1(ω) and ε2(ω), and absorption coefficient α(ω) for CuIn1-xAlxSe2. Rybicki and Lightman 41 state that the absorption coefficient indicates the travelling displacement a photon could have before being absorbed by the materials. It can be obtained directly from the dielectric function. As the experimental gaps of the CuIn1-xAlxSe2 when x is in between 0 and 1 are unknown, the scissor operation was not applied for these compounds in considering the consistency.

The evolutions of optical responses at perpendicular and parallel directions conform to each other. The decreasing imaginary part ε2(ω) as well as the right shift of the peaks as x increases from 0 to 1 when the incident energy is lower than 2.5 eV can be interpreted as the decreasing ability to absorb some light wavelengths. This can be correlated to the enlarged band gap as x increases. Incident photons with high enough energy can be absorbed by large band gap material to trigger the electrons transition, while they are mainly wasted in small band gap materials through heat relaxation.

The static dielectric function ε1(0) and frequency-dependent ε1(ω) are related to the band gap energy since they are obtained from the imaginary part ε2(ω). The Penn model 42 exhibits an inverse relation between ε1(0) and the square of the band gap energy Eg, 𝜀 + (0) ≅ 1 + (ћω/𝐸g) ! , where ћω is the plasma energy. From Fig. 6a, ε1(0) decreases from 17 to 8.5 as x increases from 0 to 1, which reflects an increase in optical band gap since higher ε1(0) gives lower band gap. The fluctuation of ε1(ω) in the energy range from Eg to around 6 eV indicates the interaction of material with the photons.

Regarding the absorption coefficient (Fig. 6b) one can note that the calculated absorption edge values are consistent with variation of band gaps at the center of the first Brillouin zone. With increasing incident photon energy, the absorption coefficient increases to reach a maximum value around 9 eV with various small peaks in between, corresponding to electron transitions from the valence band to the conduction one. In the energy range of [0;6] eV, which b a corresponds to the most intense part of sun radiation, the absorption coefficient decreases as x increases. The absorption peaks experience a right shift in the range [0;4] eV and a left shift in the range [4;6] eV. Above 6 eV, the absorption coefficients of all the compounds increases sharply to their maximum; CuAlSe2 absorption coefficient is the highest. Spectroscopic limited maximum efficiency (SLME) has been recognized as a useful method to characterize the PV abilities of materials 43 . For a real solar cell, the theoretical maximum efficiency depends on the thickness of the absorber layer 44 . In the energetic domain [1.5;3.1] eV, the maximum efficiency η of CuIn1-xAlxSe2 has been calculated with respect to the absorber thickness, as shown in Fig. 6c. The substituted system displays a decreasing η as x goes from 0 to 1 at any given thickness. The maximum one corresponds to pure CuInSe2, whereas among all the CuIn1-xAlxSe2 alloys, that with x equals 0.25 has the highest conversion efficiency. Accordingly, CuIn0.75Al0.25Se2 is a promising PV material for the fabrication of solar cells in considering the good trade-off between a large enough band gap and a relatively high conversion efficiency.

Topological properties

To unravel the possible relationship between optical properties and topological properties of the electron density, the influence of aluminum substitution has also been investigated on topological properties. Fig. 7 shows the electron density and Laplacian distributions for the optimal indium substituted compound, CuIn0.75Al0.25Se2, in both planes (110) and (-110). Electron density diagrams serve as a complementary tool to achieve a proper understanding of the electronic structure. Bader's theory defines critical points that are saddle points in the charge density on the zero-flux surfaces. Bond critical points are critical points on bond paths, the paths between bonded atoms along which the charge density is at a maximum with respect to that along any neighboring paths. It is found that there is a correlation between the values of the electron density at the bond critical points (BCPs) and the interaction strength between atoms 45 . The properties at the BCPs between Cu, In, Al nuclei and Se ones (M-Se) constitute our main subject of interest.

Besides, the Laplacian of the electron density ( ∇ ! 𝜌 :;< ) at BCPs provides details analogous to those found in relation to ρBCP. According to the local virial theorem, the Laplacian can be expressed in function of the kinetic GBCP and potential VBCP energy density as ћ ! 4𝑚 ⁄ ∇ ! 𝜌 :;< = 2𝐺 :;< + 𝑉 :;< 46,47 . The ionic character of a material can be related to charge transfer between cations and anions while covalent character is related to charge sharing. In the present case, there are Cu-Se, Al-Se and In-Se interactions leading to the building of crystal orbitals. It is hence essential to find out their featuring attributes. defined along the "n-b-n" path (where "n" are the nuclei and "b" is the bond critical point), the electron density ρBCP, the Laplacian ∇ ! 𝜌 :;< , local electron kinetic GBCP, potential VBCP and total HBCP energy densities, as well as the parameters related to the bonding characters at the nonequivalent bond critical points in unit cell. The results presented below concern the two tetrahedron sub-structures, structure I and structure II (inset of Fig. 9), which build up all the compounds. Structure I is composed of four types of atoms, namely Cu, In, Al and Se, while structure II is composed of three types of atoms, namely Cu, In, Se or Cu, Al, Se. For x equals to 0 and 1 only structure II is present whereas for x equals 0.5 only structure I is present. Structure I and structure II both exist when x is 0.25 and 0.75. we use the asterisk sign (*) to represent the longest bond length for each pairwise M-Se. Fig. 8 depicts the bond lengths and bond angles for each pairwise M-Se, both being determined from "n-b-n" path, where "n" are the nuclei and "b" is the bond critical point. As it shows, for a given value of x the bond lengths and angles behave oppositely. For example, M-Se* has a longer bond length and a lower angle than M-Se. In-Se has the smallest angle which accords with its longest bond length in comparison to Cu-Se and Al-Se. This dependence on bond length and angle provides a gauge to evaluate the bond interaction in Cu-chalcogenide materials.

Table 3

Topological and energetic properties of each pairwise M-Se. Associated units: bond length (Ǻ) defined along the "n-b-n" path (where "n" are the nuclei and "b" is the bond critical point), ρBCP (10 -2 e•Å -3 ), 𝛁 𝟐 𝝆 𝐁𝐂𝐏 (10 -2 e•Å -5 ), GBCP, VBCP, HBCP (10 -2 kJ/mol per atomic unit volume), and HBCP/ρBCP (kJ/mol per electron). I and II correspond to two different tetrahedron structures (see text), as depicted in the inset of The sign of ∇ ! 𝜌 :;< is a useful tool to distinguish the types of interactions between atoms. From a topological point of view, pure closed-shell (CS) charge-depletion interactions at BCP are marked with a low electron density and a positive Laplacian in ionic crystals (e.g. ionic bond, hydrogen bond and van der Waals interactions).The BCPs are characterized by the relatively large value of the kinetic energy density and a positive value of the total energy density 48 . By contrast, pure shared-shell (SS) exhibits a negative Laplacian with locally b a concentrated electron density (e.g. covalent and polar bonds), dominated by the potential energy density, the total energy being negative. These two types of definition are not enough though to describe the bonding properties since most crystals do not belong to these extreme prototypes. Espinosa et al. 21 divided the bond types into three classes by using the adimensional |VBCP|/GBCP ratio and introducing the bond degree (BD=HBCP/ρBCP). In addition to the pure CS (|VBCP|/GBCP<1) and SS (|VBCP|/GBCP>2), another region named transit CS with positive Laplacian and negative total energy density (1<|VBCP|/GBCP<2) has been defined, as plotted in Fig. 9, that shows the BD versus |VBCP|/GBCP. The interactions of Cu-Se, In-Se, and Al-Se all lie in the transit CS zone, between the typical ionic and covalent bonds. Our results are in agreement with Saeed et al. 22 ones that render the covalent character with a weak polarization between Al, Ga, In cations with S, Se and Te anions.

Since for the tetrahedron structures I and II, electrostatic interactions between metal and Se ions are nonequivalent due to the aluminum substitution, the corresponding bonds have been divided into M-Se and M-Se* according to their bond lengths, as described above. In structure I, as x increases, the absolute magnitudes of the bond degree, i.e. covalence degree (CD) when HBCP<0 21 , of Cu-Se and Cu-Se* first decrease (from x=0 to x=0.5), then increase to a maximum (at x=0.75) and finally decrease as x goes to 1. While in structure II, the CD of Cu-Se and Cu-Se* first increase to a maximum (at x=0.25), then decrease sharply to a minimum (at x=0.5) and goes up again to x=1. The CD maxima for Cu-Se/Se* in structure I and II occur at 0.75 and 0.25, respectively, as well as in In-Se case. As x increases, the CD of In-Se increases while that of Al-Se decreases, and In-Se has stronger bond interaction compare with that of Cu-Se/Se* and Al-Se/Se* for each x value due to its higher covalent degree. This is consistent with Wada's results 49 . In addition, it is noticeable that the evolution between covalent degree (i.e. bond degree) and |VBCP|/GBCP ratio agrees with each other, i.e. the higher the degree, the bigger the ratio. Negative total energy density is the sign of nuclei attraction. When two atoms share electrons to form a covalent bond, bond degree can be interpreted as the intensity of bond attraction. The |VBCP|/GBCP left shift of Al-Se and |VBCP|/GBCP right shift of In-Se as x increases reflect to a stronger nuclei attraction of In-Se than Al-Se. As listed in Table 3, the accumulation of electron density at BCP along each bond path increases almost linearly as the bond length decreases. Small bond length means large overlapping space in electron cloud. Generally, any mechanism that increases electron density in material also increases refractive index. The ease of electron transition from Se 2-anion full valence band to empty metal conduction band determines the optical response in the presence of incident photons. In-Se, with its largest bond degree compared to that of Cu-Se/Se* and Al-Se, contributes the most to the high absorption coefficient of CuInSe2 since the Se 4p and In 5s+5p orbitals are prominent in the valence and conduction states, respectively, near the Fermi level. These orbitals contributions may also explain the much lower band gap of CuInSe2 than that of CuAlSe2. In effect with increasing aluminum content in CuIn1-xAlxSe2, the In 5s+5p orbitals contribute less and less to the conduction band, while the contribution of Al 3s+3p is gradually increasing. Moreover, since the Al orbitals contribution to the DOS increases, the decreasing covalent bond degree of Al-Se prevails over the In-Se increasing one, which leads to an overall decrease of the covalent bond degree. Interestingly, a correlation can be made between the maximum covalent bond degrees in Fig. 9 and the partial density of states in conduction band in Fig. 4, that gives a basis to evaluate the optical response. In Fig. 4c (x=0.25), the peaks at 1.6 eV and 3.2 eV of the conduction band, which are mainly contributed from Se 4p and In 5s orbitals, and In 5p and Cu 4s ones respectively, can be correlated to the covalent bond degree maxima of Cu-Se/Se* and In-Se in structure II (see Fig. 9), originating from the transitions of electrons between Se 4p and In 5s+5p/Cu 4s hybrid orbitals. Likewise, the maxima of Cu-Se/Se* and In-Se in structure I can be correlated to the peaks in Fig. 4e (x=0.75). These results show that, pairwise In-Se contributes the most to electron transition compare to other pairwises for the In-substituted compound.

Investigation of CuIn0.75Al0.25Se2 under strains

From the above results, CuIn0.75Al0.25Se2 compound has been chosen as target to perform the biaxial strain study. The evolution of atomic geometry parameters under compressive (-) and tensile (+) in-plane strain is shown in Fig. 10. As expected, due to the constant volume lattice relaxation, there is a linear decrease of c/c0 and c/a in the strain range [-4%;+4%]. In the meantime, u(Se), which is the free x-coordinate of Se in the structure decreases linearly as can be seen from the fitting curve in Fig. 10c.

Optical properties

Band gap, absorption coefficient and dielectric functions along perpendicular (xx) and parallel (zz) directions to c-axis of CuIn0.75Al0.25Se2 under strain are shown in Fig. 11 and Fig. 12. Band gap decreases for both compressive and tensile strains, and it decreases faster under compressive strains than under tensile ones. In-plane (αxx) and out-of-plane (αzz) absorption coefficients behave oppositely under strains. No evident dispersion of α(ω) occurs when incident energy is below 1.6 eV, except a slight decrease for αxx and a slight increase for αzz as c c a a b strain goes from compressive to tensile. For incident energy above 1.6 eV, strains favour the dispersions with significant rise of αxx and fall of αzz from compression to tension. As it shows, positive strain, which lengthens the lattice a and shortens the c one, provides a bigger chance to absorb incident photons along the xx direction than along the zz one. Conversely, the negative strain will result in a bigger chance to absorb incident photons in the zz direction than in the xx one. This is consistent with the essay to extend the optical path for increasing the effective absorption in photon management as explained by Wang et al 50 .

As expected, since absorption coefficient and dielectric functions are directly related quantities, similar behavior is observed for dielectric functions in Fig. 12. Inverse behaviours of the dielectric functions along the xx and zz directions as biaxial strains are applied have been observed. Imaginary part of the dielectric function depicts the response of the compound to an outer electromagnetic field. Increased ε2_xx indicates a strengthened inner electric field, which is formed by the free minority carriers; on the contrary, declined ε2_zz shows a weakened electric field. Considering the DOS (Fig. 4c), the two major peaks observed for the imaginary part at 1.6 eV and 3.2 eV (Fig. 12a) can be attributed to the electron transition from Se 4p to In 5s and In 5p/Cu 4s+3d, respectively. These two peaks are blue-and red-shifted for the in-plane ε2_xx component while red-and blue-shifted for the out-of-plane ε2_zz one as strain goes from negative to positive values. Real parts of dielectric functions (ε1_xx and ε1_zz) show the same trend as that of the imaginary parts. Static dielectric constants ε(0) shows larger deviation under strain for out-of-plane component than for the in-plane one.

Topological properties

To better understand the evolution of electronic and optical properties as biaxial strain is applied, the topological properties at bond critic points have been studied, as presented in Fig. 13. Both nonsymmetrical structure I (Fig. 13a) and structure II (Fig. 13d) have been considered separately. As it shows, strain has a weak influence on the atomic interactions. All pairwises interactions lie in the transit closed-shell zone, i.e. 1<|VBCP|/GBCP<2, and thus correspond to covalent bond. Narrowed and broadened dispersion of covalent degrees and |VBCP|/GBCP ratios, which can be associated to "rigid" and "loose" topological cells, respectively, are observed under compressive strains for the former and tensile strains for the latter. These two different topological cells explain the anisotropic response to electromagnetic field, such as dielectric functions. Scattered bond degree distribution reflects a large extent in bond polarity. This distribution provides information on the preferable site of maximal shared electron density and on the ease of electron transition. As we know, negative total energy density is a sign of nuclei attraction. The bare change in |VBCP|/GBCP ratios with negligible differences among Cu-Se/Se* and In/Al-Se imply a uniform and stable bonding characteristics at compressive strain, which echoes to "rigid" cell. Similarly, the decrease in |VBCP|/GBCP ratios with considerable differences assigns large deviation in bonding uniformity and intensity at tensile strain, which is in line with "loose" cell. One should underline that, the upward and downward changes of the bond degrees of Cu-Se/Se*and In-Se/Se* in structure II when strain is below -2% behave discordantly to those of the other strain ranges. This abnormal evolution of the bond degrees, which is also noticeable for |VBCP|/GBCP ratios, suggests the importance in strain-controlling during PV cell manufacturing.

As bond lengths and bond angles in Fig. 13c and 13f show, strain modifies the bond configuration by adjusting the displacement of bond critical point, as well as the distance from terminal atoms (M, Se) to bond critical point. In each crystal cell, there are four-centered structures of type I sandwiched along the c-axis by two structures of type II at each side. Structural parameters depart from the free-strained ones when unit cell is subjected to biaxial strain. Among these deviations, those pertaining to long bond lengths are higher than those pertaining to the short ones. For example, the In-Se bond length has varied more significantly than the Al-Se one, which is barely changed throughout. As we noticed previously, bond length evolves oppositely to bond angle. As to the model "nM-b-nSe", small bond angle means shorter distance from BCP to M cation than to Se anion. Displacement of BCP originates from the nuclei attraction will leads to a bond path deflection. High deflection shows large overlapping area in between shared electrons. BCP close to a M cation indicates large overlapping region in electron density with a strong electronic attraction from the M cation. Electron transition from Se anion to M cation becomes easier. By contrast, BCP close to a Se anion indicates small overlapping region in electron density, thus it is more difficult for electron to transit from Se anion to M cation. Bond angle can then be a sign of electron transition feasibility. The downward change of the band gap under both compressive and tensile strains (Fig. 11a) can be paralleled with the decrease of the bond angles; the faster decrease of the band gap at compressive strain is concomitant with its faster reduction of the bond angle than at tensile strain, particularly for pairwise In-Se. As seen above, for CuIn0.75Al0.25Se2, i) pairwise In-Se in structure II are mainly responsible for the electron transitions when incident energy is above 1.6 eV, ii) the peaks of the in-plane dielectric function are attributed to the electron transition between Se 4px/y-Cu 3dxy and Se 4px/y-Cu 3dx2-y2, while those of the out-of-plane dielectric function are assigned to the electron transition between Se 4pz-In 5s and Se 4pz-In 5pz, iii) under both tensile and compressive strains, covalent degrees of In-Se/Se* decrease. As we know, covalent degree stands for atomic interaction intensity and strong interaction leads to high absorption coefficient. αxx goes down under compressive strain (lattice a decreases), whereas αzz behaves oppositely. The decrease of αxx and αzz under compressive and tensile strain, respectively, can be correlated to the negative shift of the covalent bond degree of In-Se. Therefore, strengthened absorption occurs along the lengthened direction, i.e. αzz at compressive strain and αxx at tensile strain, which could be due to the long transportation path for incident photon providing a big chance to be absorbed before it escapes. This is consistent with the downward evolution of bond angle under strains, which gives high feasibility for electron transition.

Conclusions

In this study, the band structure, density of states, absorption coefficient, and dielectric functions of CuIn1-xAlxSe2 have been calculated by using the full potential linear augmented plane wave (FP-LAPW) method. CuIn1-xAlxSe2 has a direct band gap, its increasing value as aluminum substitutes indium contrasts with their continuously declining absorption coefficients. CuIn0.75Al0.25Se2 has been recognized as the optimal alloy in terms of the compromise between band gap and absorption coefficient. Based on electron density ρBCP, and the |GBCP|/VBCP, HBCP/ρBCP ratios at BCP, the nonequivalent M-Se bond interactions all lie in the transit closedshell zone and show intermediate polar covalent interactions with positive Laplacian and negative total energy (1<|VBCP|/GBCP<2). The In-Se bond bears the largest covalent degree. As x increases, the covalent degree of pairwise Al-Se decreases. The connection established between bond length, bond angle and bond degree for each pairwise M-Se, can be used as a tool to predict the ease of electron transition. The absorption coefficients in perpendicular and parallel direction to c-axis increase and decrease, respectively, as strain goes from compressive to tensile. Shrinkage and expansion of the dispersions in bond degree due to compressive and tensile strains, respectively, are responsible for the anisotropic optical responses. Inverse evolution for perpendicular and parallel components of dielectric function and absorption coefficient under strain originates from the bond asymmetry. The decrease of band gap with respect to strain is asymmetric, it decreases faster under compressive strain than tensile one. The whole set of results obtained in this work show the strong anisotropic impact of the biaxial strains on the topological and optical properties of CuIn1-xAlxSe2. This suggests managing the strains as best as possible during the PV cell manufacturing, especially during epitaxial growth.

Fig. 1

 1 Fig. 1 Structure optimization of CuInSe2 (a) (b) and CuAlSe2 (c) (d).

Fig. 2

 2 Fig. 2 Lattice parameters and their fitting results as aluminum content x varies.

  The band structures of CuIn1-xAlxSe2 alloys (x=0, 0.25, 0.5, 0.75, and 1) are plotted in Fig.3along the high symmetry directions in k-space. The dashed horizontal line in band structure shows the position of Fermi level energy (EF), and was set to zero as reference in this work. The values of band gap Eg can be obtained from the energy difference between conduction band minimum (CBM) and valence band maximum (VBM). The highest valence band states originate from a threefold-degenerate p orbital. Both the CBM and the VBM occur at Γ point, the center of the first Brillouin zone, indicating a direct band gap. As shown in Fig.

Fig. 3

 3 Fig. 3 Band structures of CuIn1-xAlxSe2 when x equals to 0 (a), 0.25 (b), 0.5 (c), 0.75 (d), 1 (e), and their corresponding band gap (f).

  2 eV, in addition to Cu d and Se p orbitals, there is a small contribution of In p and Al p orbitals in the mixing. In the region of [-5.4;-4.6] eV, we observe a comparatively strong hybridization of Cu s and Se p states. Conduction band between 1.5 eV and 2.5 eV is

Fig. 4

 4 Fig. 4 Total density of states (a) and partial density of states of CuIn1-xAlxSe2when x equals to 0 (b), 0.25 (c), 0.5 (d), 0.75 (e), and 1 (f).

Fig. 5

 5 Fig. 5 Calculated (solid lines) and experimental (dash lines) results of dielectric functions of CuInSe2 (a) and CuAlSe2 (b) in the perpendicular (xx, black lines) and parallel (zz, red lines) directions to c-axis. The experimental data are taken from Ref. 39,40 .

Fig. 6

 6 Fig. 6 Dielectric functions (a) and absorption coefficient (b) for CuIn1-xAlxSe2 at different x values. The solid lines and dashed lines represent the perpendicular (xx) and parallel (zz) directions to c-axis, respectively. Calculated maximum efficiency (c) with respect to the film thickness in the energetic domain [1.5;3.1] eV.

Fig. 7

 7 Fig. 7 Electron density (top panels) and Laplacian (bottom panels) distribution of CuIn0.75Al0.25Se2 in plane (110) (a) (c) and plane (-110) (b) (d). b, c, and r mean bond, cage, and ring critical points. Green lines represent the electron density ρ, magenta and blue lines represent the negative and positive ∇ ! 𝜌, respectively.

Fig. 9 .

 9 Fig.9. The asterisk sign (*) means the longest bond length for each pairwise M-Se.

Fig. 8 9 .

 89 Fig. 8 Calculated bond lengths (a) and bond angles (b) in CuIn0.75Al0.25Se2. I and II correspond to two different tetrahedron structures (see text), as depicted in the inset of Fig. 9. The asterisk sign (*) stands for the longest bond length for each pairwise M-Se. Both the bond lengths and angles are determined along the "n-b-n" path, where "n" are the nuclei and "b" is the bond critical point, respectively.

Fig. 9

 9 Fig. 9 Values of bond degree (BD=HBCP/ρBCP) vs. |VBCP|/GBCP at BCPs of M-Se as x varies from 0 to 1. Structures I (Cu/In/Al-Se, a) and II (Cu/In-Se or Cu/Al-Se, b) represent two different constitutive parts of every unit cell. The asterisk sign (*) stands for the longest bond length for each pairwise M-Se.

Fig. 10

 10 Fig. 10 Geometric characteristics of CuIn0.75Al0.25Se2 and their evolution under biaxial strains (a-a0)/a0. Unit cell (a), green, red, yellow and magenta atoms represent In, Cu, Se and Al, respectively. lattice constant ratio c/c0 and c/a (b) and internal parameters u(Se) (c).

Fig. 11

 11 Fig. 11 Band gap vs. biaxial strain(a-a0)/a0(%) (a) and absorption coefficient in perpendicular (xx) and parallel (zz) directions to c-axis vs. incident energy for various biaxial strains (a-a0)/a0(%) (b).

Fig. 12

 12 Fig. 12 Imaginary part (a) and real part (b) of dielectric functions at perpendicular (xx) and parallel (zz) directions to c-axis vs. incident energy for various biaxial strains (a-a0)/a0(%).

Fig. 13

 13 Fig. 13 Topological properties at bond critic points of structure I (a) and structure II (d) vs. biaxial strain (a-a0)/a0(%). (b) and (e): ratio of potential energy density to kinetic energy density |VCP|/GCP and bond degree BD; (c) and (f): bond length BL and bond angle.

Table 1

 1 Lattice parameters a and c and bond lengths d(Ǻ) of CuInSe2 and CuAlSe2. u(Se) is the free x-coordinate of Se in the chalcopyrite unit cell.

			GGA functionals	a	c	u(Se)	c/a	dCu-Se	dIn-Se	dAl-Se
		Exp. 32		5.7815 11.6221 0.2426 2.0102 2.5326 2.4833	-
	CuInSe2	Calc.	PBE PBEsol	5.8798 11.7672 0.2688 2.0013 2.5723 2.5236 6.0369 11.9387 0.2814 1.9776 2.4430 2.6551	--
			WC	5.8775 11.7565 0.2683 2.0002 2.4845 2.6089	-
		Exp. 33		5.6060 10.9010 0.2600 1.9443 2.4383	-	2.3729
	CuAlSe2	Calc.	PBE PBEsol	5.5307 10.7536 0.2506 1.9443 2.3748 5.5983 11.0090 0.2571 1.9443 2.3878	--	2.3709 2.4340
			WC	5.6194 10.9261 0.2506 1.9443 2.4129	-	2.4089

Table 3

 3 lists the bond distances

	a	b	c

Acknowledgements

This work is financially supported by China Scholarship Council (CSC). This work was granted access to the HPC resources of the Centre Informatique National de l'Enseignement Supérieur (CINES), Montpellier, France under allocation A0050806881 made by the Grand Equipement National de Calcul Intensif (GENCI). It was also granted access to the HPC resources of Aix-Marseille Université financed by the project Equip@Meso (ANR-10-EQPX-29-01) of the program "Investissements d'Avenir" supervised by the Agence Nationale de la Recherche.

Conflicts of interest

There are no conflicts of interest to declare.