Modeling of the Interminiband Absorption Coefficient in InGaN Quantum Dot Superlattices

Abstract: In this paper, a model to estimate minibands and theinterminiband absorption coefficient for a wurtzite (WZ) indium gallium nitride (InGaN) self-assembled quantum dot superlattice (QDSL) is developed. It considers a simplified cuboid shape for quantum dots (QDs). The semi-analytical investigation starts from evaluation through the three-dimensional (3D) finite element method (FEM) simulations of crystal mechanical deformation derived from heterostructure lattice mismatch under spontaneous and piezoelectric polarization effects. From these results, mean values in QDs and barrier regions of charge carriers’ electric potentials and effective masses for the conduction band (CB) and three valence sub-bands for each direction are evaluated. For the minibands’ investigation, the single-particle time-independent Schrödinger equation in effective mass approximation is decoupled in three directions and resolved using the one-dimensional (1D) Kronig–Penney model. The built-in electric field is also considered along the polar axis direction, obtaining Wannier–Stark ladders. Then, theinterminiband absorption coefficient in thermal equilibrium for transverse electric (TE) and magnetic (TM) incident light polarization is calculated using Fermi’s golden rule implementation based on a numerical integration into the first Brillouin zone. For more detailed results, an absorption coefficient component related to superlattice free excitons is also introduced. Finally, some simulation results, observations and comments are given.


Introduction
In photonics and optoelectronics, many devices are based on photon absorption-e.g., photodetectors, solar cells, and so on-or are influenced in a weak or strong manner by it-e.g., optical waveguides, fibers, modulators, and so on.This very important material capability strongly depends on the semiconductor crystal type used, and it is tightly linked to the specific device component structure, as well.Indeed, its values and features change for cubic zincblende (ZB), e.g., gallium arsenide (GaAs) or indium arsenide (InAs), with respect to hexagonal wurtzite (WZ) semiconductors, e.g., binary III-nitride compounds, such as gallium nitride (GaN) or indium nitride (InN) or their ternary or quaternary alloys, such as InGaN.Moreover, the photon absorption, as all semiconductor optical properties, changes as a function of the considered structure, either the bulk or nanostructure-e.g., quantum well (QW), multi-quantum well (MQW), QW superlattice (QWSL), nanowire, single QD, multi-QD (MQD) or QDSL in nanophotonics and nano-optoelectronics.
In this paper, our interest has been focused on self-assembled QDSL and InGaN WZ materials [1][2][3][4].A semi-analytical computational approach for QDSL minibands based on the model developed by Lazarenkova and Balandin [5] and the interminiband absorption coefficient has been implemented in order to apply it to the InGaN semiconductor nanostructures.Compared toprevious works with a similar miniband calculation methodology [2][3][4][5][6][7], our method introduces, into minibands' analysis, the strain influence on charge carriers'electric potentials and effective masses evaluated from k¨p theory under spontaneous and piezoelectric polarizations [1,8], considering the minibands' modification due to an internal electric field component along the semiconductor crystallographic polar axis (F c ), as well.Furthermore, we have also analyzed the miniband formation in WZ valence sub-bands, labeled as A, B and C, from that with the upper-lower energy level.Into the QDSL interminiband absorption coefficient estimation, minibands' non-parabolicity, light polarization, as investigated in previous works for bulk [9] and QW [9,10] structures, and strain influence on the interband momentum matrix element have been included.The last issue introduced into the model has been considered in a few previous works.For example, for WZ GaN and aluminum nitride (AlN), it is evaluated through density functional theory (DFT) first-principles calculations and is interpreted in the contest of k¨p theory [11].Instead, here, it is derived using the k¨p formulation for Kane parameters-K η x , K η y and K η z for η-conduction band (CB) charge carrier transition, where η represents one of A, B or C valence sub-bands-without strain, in which the strain-independent effective masses are replaced with those modified by strain (see Appendix C).Furthermore, for more detailed results, an absorption coefficient component related to superlattice free excitons (see Appendix D) is also introduced.
A superlattice was introduced for the first time in 1970 by Esaki and Tsu with a work on a 1D periodic potential for electrons formed by a series of alternate regions of different materials (dissimilar alloy compositions or unlike impurity concentrations), with a period size of several nanometers, but shorter than the electron mean free path [12].This structure allowed energy level minibands for charge carriers (absent in bulk structures) to be present [13].Minibands' electric and optical properties can be engineered by changing the semiconductors and/or region sizes of the superlattice.Thus, this heterostructure can become useful for many applications in which specific desired energy levels or bands are needed-e.g., in full solar spectrum photovoltaic devices [14,15].A one-dimensional (QW) superlattice is the most investigated and used superlattice type, due to its simpler physical-mathematical modeling and more mature fabrication technology (generally based on epitaxial growth).Indeed, photodiodes [16], solar cells [14,15], LEDs [17][18][19], lasers [20] and other devices in which one or more regions-e.g., absorption region, reflector, gain active region, etc.-are formed by a QWSL have been fabricated or investigated.On the contrary, two-or three-dimensional (2D or 3D) superlattices' (nanowire and QD superlattices, respectively) growth technology is more complex, and the evaluation of its optical properties is generally characterized by large time and memory consumption, although a plane wave expansion is used-e.g., where each 3D charge carrier wavefunction evaluated for a given q vector of the respective 3D superlattice reciprocal space is expanded with 2197 plane waves [21].
On the other hand, in the same years, studies on InGaN growth methods-e.g.,the electron beam plasma technique [22] or metal-organic chemical vapor deposition (MOCVD) [23]-led to a progressive increase of interest on this ternary alloy as a light source and absorption material.In the 1990s, the first InGaN LEDs with different emission wavelengths were demonstrated [24][25][26].With the possibility to obtain single-crystalline InN layers on sapphire (Al 2 O 3 ) substrate by molecular beam epitaxy (MBE) [27] and to measure its low band gap [28], in the 2000s, investigations on InGaN absorption and the use of this semiconductor in solar cells increased [29] for both bulk structures [30,31] and nanostructures [15,16].Then, InGaN become a highly interesting material for high performance solar cells.Recently, researches have focused their interests on more specific InGaN features that can influence photovoltaic devices' operation [32,33].

InGaN Semiconductor Features
InGaN is a ternary alloy semiconductor derived from the two binary III-nitride compounds, GaN and InN.It is thermodynamically stable in the hexagonal WZ crystal structure.Indium (In) content x of the In x Ga 1´x N compound, variable from 0 (GaN) to 1 (InN), influences its mechanical, electric and optical properties.Numerical values of the parameters that describe all of these properties can be obtained with a linear interpolation (Vegard's law) of those of GaN and InN, except for band gap E g (E In x Ga 1´x N g ) and spontaneous polarization P (P In x Ga 1´x N ) for which their own corrective bowing parameters are needed in a parabolic interpolation.The spontaneous polarization P is a polarization effect along the polar axis crystallographic direction [0001] (i.e., the nomenclature for the WZ lattice) absent along all directions into the plane of the other two crystallographic axes (called here "in-plane"), which characterizes all WZ semiconductors.Thus, the In x Ga 1´x N band gap varies according to Equation (1) from the value for InN E InN g = 0.78 eV to that for GaN E GaN g = 3.51 eV with the bowing parameter E Bow g = 1.4 eV.Furthermore, spontaneous polarization P In x Ga 1´x N can be described by an equation similar to Equation (1) in which band gaps E InN  g , E GaN g and E Bow g are replaced with P InN , P GaN and P Bow , respectively.All physical parameter values used in this paper are taken from [34,35] and tabulated in Table 1.
Table 1.GaN and InN physical parameters taken from [34,35] and used in this paper.For In x Ga 1´x N, they are obtained as linear or parabolic (with theirown bowing parameter) interpolation of those of GaN and InN with respect to In content x.The general superscript tag "α" substitutes labels "GaN", "InN" or "Bow" (bowing) used in the text for GaN and InN material and bowing parameters, respectively.´5.5 ´5.5

Physical Parameters
Moreover, In x Ga 1´x N lattice constant values c In x Ga 1´x N (along the polar axis) and a In x Ga 1´x N are a linear function of those of involved binary compounds GaN (c GaN and a GaN ) and InN (c InN and a InN ).Lattice mismatches ε c (along the polar axis) and ε a (in-plane mismatch) of the In x Ga 1´x N/In y Ga 1´y N heterostructure are evaluated as in Equation ( 2): In an In x Ga 1´x N/In y Ga 1´y N heterostructure with polar c-plane growth, i.e., along the InGaN polar axis [0001], the in-plane lattice mismatch ε a maximum value is 11% and is obtained for In contents x = 0 and y = 1 (InN growth on the GaN layer) [36].Mismatch during hetero-epitaxial growth induces a compressive or tensile stress, and then, the epilayer grows as a pseudomorphic layer; strained films with deformation of their own lattice and energy band structure [34].Stress increases with layer thickness until a mechanical relaxation process occurs.Two types of relaxation mechanisms can occur.In the former, hetero-stress is reduced with a 3D cluster formation (Stranski-Krastanov (SK) self-assembled QD growth mode) beyond a previously-grown thin pseudomorphic layer called the wetting layer.In the latter, a misfit dislocation formation relaxes the grown layer stress.While in some applications, the former can be desired, the latter induces low quality layers with a high dislocation density, particularly on hetero-interfaces [14,15].Zhao et al. [37] and Pristovsek et al. [38] described theoretically and experimentally the InGaN growth with metal-organic vapor phase epitaxy (MOVPE) on a GaN buffer layer having a [0001] growth direction.Their attention was focused on In content and the growth process parameters' influence on the critical thicknesses for two relaxation mechanisms, in order to establish what occurs in each heterostructure.
The strain effect is present in all InGaN devices that can be fabricated only as heterostructures, because unfortunately, InGaN native homogeneous substrate fabrication is today still a challenge for semiconductor growth technology.Then, an InGaN layer is grown almost always in hetero-epitaxial mode on GaN that is characterized by different lattice constants.In turn, GaN can be also grown in hetero-epitaxial mode as a buffer layer on a different substrate type with a different lattice constant, sometimes called a "foreign substrate", such as sapphire with a 16% mismatch in crystal constant [39], silicon (Si), silicon carbide (SiC) [36] or GaAs [26,39], introducing into devices lattice mismatch stress and, consequently, mechanical strain and high dislocation defect density.Strain and, above all, dislocation defects can be reduced with some techniques based on removing the foreign substrate on which GaN was grown (sapphire substrateremoved through thelaser lift off (LLO) technique and etching by diluted hydrochloric acid (HCl) solution [40] or GaAs removed mechanically [26,39]), obtaining a GaN freestanding substrate.Recently, an ammonothermal GaN growth technique has allowed the direct fabrication of a native GaN substrate, and it has been also introduced in high quality LED commercial production [41].
GaN and InGaN are generally grown along the polar axis [0001], due to theirmature technology with respect to that of other growth directions (non-polar a-plane 1120 ( , m-plane 1100 ( or other planes normal to semipolar growth directions).InGaN can be grown directly also on Si, as already demonstrated in [42].Further, the InGaNfilms or QDs strain, as in all WZ semiconductors, induce a piezoelectric polarization, as well.Spontaneous and piezoelectric polarizations in the InGaN device regions led to a built-in electric field that modifies their operation and produces a spatial electron-hole separation with, sometimes, the consequent decrease of the absorption coefficient: the quantum confined Stark effect (QCSE) in devices based on nanostructures.
As that of its binary III-nitride compounds, for the InGaN electronic band structure at the center of the Brillouin zone (the Γ symmetry point in which all orthogonal components q x , q y and q z respectively along the x, y and z directions of the reciprocal space vector q are zero), we consider only CB and the three A, B and C valence sub-bands that are involved in electronic and optical processes as optical absorption.Their doubly degenerate (two states for different charge carrier spin Ò (up) and Ó (down)) energy levels E A , E B and E C shown in Equation (3) and their wavefunctionsu A 1{2 , u B 1{2 and u C 1{2 in Equation (4) (VB Bloch lattice functions for semiconductor under spin-orbit and crystal field splitting) are eigenvalues and eigenvectors of the diagonalized WZ one-hole Hamiltonian matrix [8,[43][44][45] without any strain effect derived by the k¨p theory of Luttinger-Kohn: In Equation ( 3), cr and so are semiconductor crystal field and spin-orbit splitting energies, respectively, and a and b are defined in Equation (5).In x Ga 1´x N cr and so splitting energies are functions of In content x and derive from Vegard's law applied on the corresponding parameters for GaN and InN semiconductors (Table 1).
Doubly-degenerate CB energy level E CB and Bloch lattice functions u CB 1{2 are shown in Equation ( 6) in which î is an imaginary unit.
In Equations ( 4) and ( 6) |Xy, |Yy, |Zy and |Sy are the three VB and CB Bloch lattice functions in absence of spin-orbit interaction.The former exhibit symmetry properties and a shape analogous to p atomic orbitals, the latter to s ones [45].

QDSL Minibands and the Interminibands Absorption Coefficient
In this section, a model is developed in order to evaluate the minibands and absorption coefficient due to electron transitions from minibands in valence Sub-bands A, B and C to those in the CB of a generic In x Ga 1´x N/In y Ga 1´y N QDSL (with x > y), as sketched in Figure 1.This is a heterostructure that could be grown on a GaN buffer on a foreign substrate (see Section 2), but in this study, the influence of the strain induced by the GaN buffer and substrate on the superlattice is neglected.This assumed condition sometimes can introduce into our results a little spectral shift with respect to the actual MQDs structure, but allows a more accurate analysis of minibands' effect on the absorption coefficient.
We assume our samples under investigation to be in thermal equilibrium.In this condition, CB and VB quasi-Fermi level energies correspond to the unique carrier Fermi level, which is assumed as placed between the minibands of the VB and CB.This represents an ideal situation in which VB states are completely full and CB states totally empty.Moreover, the Fermi level is assumed constant throughout the whole structure.Both of these assumptions are not actual in operative conditions for a QDSL used in a device for two reasons.First of all, even in the absence of any internal or external electric field, just a very low light absorption, a process that we are analyzing in this paper, induces a splitting of electron and hole quasi-Fermi energies (the condition of quasi-equilibrium).Secondly, QDSL is generally inserted between other structures-e.g., between a region doped with acceptor impurity (p-type) and another with donor impurity (n-type) of the p-n junction-that spatially modifies the CB and VB edge with respect to the Fermi energy level.In both cases, the occupation of the charge carriers' states will change, so reducing the light absorption.However, the QDSL absorption spectrum in thermal equilibrium gives a theoretical indication, and, e.g., its values could be used in a model for optical absorption based on minibands' capture and emission rates [46], which is already considered the charge carrier states' occupation.
Photonics 2016, 3, 5 6 of 25 splitting of electron and hole quasi-Fermi energies (the condition of quasi-equilibrium).Secondly, QDSL is generally inserted between other structures-e.g., between a region doped with acceptor impurity (p-type) and another with donor impurity (n-type) of the p-n junction-that spatially modifies the CB and VB edge with respect to the Fermi energy level.In both cases, the occupation of the charge carriers' states will change, so reducing the light absorption.However, the QDSL absorption spectrum in thermal equilibrium gives a theoretical indication, and, e.g., its values could be used in a model for optical absorption based on minibands' capture and emission rates [46], which is already considered the charge carrier states' occupation.Another approximation used here is on QDs' geometrical structure.Indeed, actual and most investigated self-assembled QDs' shape is sometimes like a hexagonal pyramid [47,48], a truncated hexagonal pyramid [1,49,50] or a lens [51], but since minibands' properties derive from a very large QD number, in the first approximation, it is possible to assume that QDSL miniband formation and properties are less influenced by the exact QDs shape when compared to their mean sizes and inter-distances.Indeed, our numerical investigations performed in isolated In0.4Ga0.6N/GaNcuboid QD with sizes 6×6×3 nm 3 reveals that the overlap value of the envelope functions is 0.77, while in [47], for the In0.4Ga0.6N/GaNQD with a truncated hexagonal pyramid shape, a value of about 0.72 has been obtained.As a result, it is possible to assume that the error on the absorption coefficient values is less than almost 7% if the truncated hexagonal pyramid QDSL is approximated with cuboid QDSL.In this context, an orthorhombic symmetry and square base is assumed for the QDs' shape (square cuboid QDs), which can simplify the superlattice modeling with the QDs' distribution regularity in the three directions [1][2][3][4][5][6][7]18,52].Indeed, with these assumptions, QDSL can be considered as a semiconductor crystal with a tetragonal lattice.Furthermore, the wetting layer influence in this work has been neglected.
It is worth outlining that the procedure adopted in this work is based on a multiphysics approach, i.e., the 3D-FEM mechanical module has been used together with the 3D-FEM piezoelectric module in order to carry out the stress distribution, and the electric field is generated by the spontaneous polarization and piezoelectric effect.Once the stress distribution and the electric field have been calculated, the energy levels in CB and in all valence sub-bands with their relative effective masses have been evaluated by applying the k•p theory (see Appendix B).At this step, the QDSL miniband formation has been carried out by applying the Kronig-Penney model to 1D periodic structures obtained by decoupling the 3D problem (see Appendix A).Moreover, the influence of built-in electric field dominant component has been evaluated as detailed in Appendix A. Finally, the absorption coefficient in thermal equilibrium (see Appendix C) has been calculated by considering the non-parabolicity of minibands, the light polarization influence and the strain dependence, as well asintroducing the absorption coefficient component related to free Another approximation used here is on QDs' geometrical structure.Indeed, actual and most investigated self-assembled QDs' shape is sometimes like a hexagonal pyramid [47,48], a truncated hexagonal pyramid [1,49,50] or a lens [51], but since minibands' properties derive from a very large QD number, in the first approximation, it is possible to assume that QDSL miniband formation and properties are less influenced by the exact QDs shape when compared to their mean sizes and inter-distances.Indeed, our numerical investigations performed in isolated In 0.4 Ga 0.6 N/GaN cuboid QD with sizes 6ˆ6ˆ3 nm 3 reveals that the overlap value of the envelope functions is 0.77, while in [47], for the In 0.4 Ga 0.6 N/GaN QD with a truncated hexagonal pyramid shape, a value of about 0.72 has been obtained.As a result, it is possible to assume that the error on the absorption coefficient values is less than almost 7% if the truncated hexagonal pyramid QDSL is approximated with cuboid QDSL.In this context, an orthorhombic symmetry and square base is assumed for the QDs' shape (square cuboid QDs), which can simplify the superlattice modeling with the QDs' distribution regularity in the three directions [1][2][3][4][5][6][7]18,52].Indeed, with these assumptions, QDSL can be considered as a semiconductor crystal with a tetragonal lattice.Furthermore, the wetting layer influence in this work has been neglected.
It is worth outlining that the procedure adopted in this work is based on a multiphysics approach, i.e., the 3D-FEM mechanical module has been used together with the 3D-FEM piezoelectric module in order to carry out the stress distribution, and the electric field is generated by the spontaneous polarization and piezoelectric effect.Once the stress distribution and the electric field have been calculated, the energy levels in CB and in all valence sub-bands with their relative effective masses have been evaluated by applying the k¨p theory (see Appendix B).At this step, the QDSL miniband formation has been carried out by applying the Kronig-Penney model to 1D periodic structures obtained by decoupling the 3D problem (see Appendix A).Moreover, the influence of built-in electric field dominant component F z has been evaluated as detailed in Appendix A. Finally, the absorption coefficient in thermal equilibrium (see Appendix C) has been calculated by considering the non-parabolicity of minibands, the light polarization influence and the strain dependence, as well asintroducing the absorption coefficient component related to free excitons of the superlattice (see Appendix D).The procedure sequence is shown in the following flow chart (Scheme 1).

Simulations, Results and Discussion
First of all, in this section, our developed model has been applied to investigate an In0.4Ga0.6N/GaNcuboid QDSL.In such a structure, dots and barrier regions band gaps are equal to 2.08 eV and 3.51 eV, respectively, while neglecting strain, spontaneous and piezoelectric polarizations.This structure will be characterized by CB (VB) miniband formation with energy levels between those of semiconductor CB (VB) edges in QDs and barrier regions.
Then, it is possible to predict spectrum values for the interminiband transition absorption coefficient included into 2.08-3.51eV (a light wavelength of 353-596 nm).We assume dot base width = 4-6 nm, dot height = 3 nm and inter-dot distances = 5 nm and = 2.5 nm.In Figure 2, we show the strain distribution of the dominant x component as evaluated by FEM simulations (see Appendix B).

Simulations, Results and Discussion
First of all, in this section, our developed model has been applied to investigate an In 0.4 Ga 0.6 N/GaN cuboid QDSL.In such a structure, dots and barrier regions band gaps are equal to 2.08 eV and 3.51 eV, respectively, while neglecting strain, spontaneous and piezoelectric polarizations.This structure will be characterized by CB (VB) miniband formation with energy levels between those of semiconductor CB (VB) edges in QDs and barrier regions.
Then, it is possible to predict spectrum values for the interminiband transition absorption coefficient included into 2.08-3.51eV (a light wavelength of 353-596 nm).We assume dot base width L β,dot = 4-6 nm, dot height L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,bar = 2.5 nm.In Figure 2, we show the strain distribution of the dominant x component as evaluated by FEM simulations (see Appendix B).

Simulations, Results and Discussion
First of all, in this section, our developed model has been applied to investigate an In0.4Ga0.6N/GaNcuboid QDSL.In such a structure, dots and barrier regions band gaps are equal to 2.08 eV and 3.51 eV, respectively, while neglecting strain, spontaneous and piezoelectric polarizations.This structure will be characterized by CB (VB) miniband formation with energy levels between those of semiconductor CB (VB) edges in QDs and barrier regions.
Then, it is possible to predict spectrum values for the interminiband transition absorption coefficient included into 2.08-3.51eV (a light wavelength of 353-596 nm).We assume dot base width = 4-6 nm, dot height = 3 nm and inter-dot distances = 5 nm and = 2.5 nm.In Figure 2, we show the strain distribution of the dominant x component as evaluated by FEM simulations (see Appendix B).After FEM 3D piezo-electro-mechanical simulations and application of the 1D Kronig-Penney model for each direction, barrier potentials for charge carriers in QDs and allowed energy minibands are found as indicated in Table 2.
Table 2.In 0.4 Ga 0.6 N/GaN QDSL electronic features for L z,dot = 3 nm, L β,bar = 5 nm, L z,bar = 2.5 nm and L β,dot = 4-6 nm without miniband splitting due to the built-in electric field.As shown in Table 2, an increase of dots' width generally reduces the QD barrier potential for charge carriers in CB and A and B valence sub-bands, while for those in C valence sub-bands, there is an opposite effect due to the strain, as can be seen from Equation (B5).This particular effect, as verified from FEM 3D piezo-electro-mechanical simulations, is the result of two contrary strain contributions in QDs and barrier regions.Indeed, while into larger dots, compressive strain and consequently an increase of the potential edge arereduced, tensile strain in barriers is applied to a more extended region, so increasing its total effect and lowering the barrier edge energy level.Between these two effects, those in barrier regions are dominant and reduce the QD barrier potential.On the other hand, this effect does not lead to a reduction of the numbers of minibands, because it is above all related to barrier regions' relative dimensions with respect to those of dots for each direction.Then, for smaller inter-dot distances with respect to dot sizes, there is a tendency to have more minibands, e.g., an increase of the number of 3D minibands can be obtained with increasing the dot size, without changing that of barriers.
In Figure 3, electron miniband dispersion in q-space for the third structure of Table 2 is shown along the two representative [[100]] (Figure 3a) and [[001]] (Figure 3b) directions (here, a double square bracket for differently labeling superlattice reciprocal space directions is used with respect to the single-square bracket for the directions of the crystal reciprocal space).Fourteen very thin minibands-energy widths with values around that of 1-12 meV-stay under the barrier continuum energy edge (red line) of both directions.Five couples of these minibands are degenerated in the Γ point and are symmetric over the plane in the Γ-M direction in the 3D reciprocal space.This is due to the tetragonal crystal structure of the square cuboid QDSLs analyzed.One of each couple of these minibands is indicated with round brackets in Table 2.
In the results shown in Figure 3, in order to evidence the features of the miniband formation, we have considered the influences of the built-in electric field in the stress distribution and then in the average values for the barrier potential and effective masses evaluated (direct influence in Equations (B1)-(B5) and indirect in Equations (B7)-(B8) of Appendix B), neglecting the direct and strong change of the Schrödinger equation (Equation (A1) of Appendix A) and charge carriers' steady states that correspond to neglecting the energy band bending typically present in WZ materials.
For the nanostructures under investigation in this work, the value of the z-component F z obtained from 3D FEM simulations is dominant with respect to F x and F y .For this reason and for a more accurate analysis of this effect, we consider the electric field only along the polar axis.Furthermore, the values of electric field F z along the central axis of cuboid QD in the z direction obtained are different for dot and barrier regions, but almost constant.For example, from FEM simulation in the second sample of Table 2, two different values of the built-in electric field, F z,dot = 0.99 MV/cm and F z,bar = ´1.02MV/cm, are achieved.Since QDs' semiconductor lattice constant is larger than that in the barrier regions, a compressive stress in QDs and a tensile one inbarriers is induced in QDSL (Figure 2) with a consequently opposite sign for F z,dot and F z,bar along the z direction.Moreover, in the numerical evaluations of F z,dot and F z,bar , the spontaneous polarization provides a very low contribution with respect to the piezoelectricity.
Photonics 2016, 3, 5 9 of 25 obtained are different for dot and barrier regions, but almost constant.For example, from FEM simulation in the second sample of Table 2, two different values of the built-in electric field, = 0.99 MV/cm and = −1.02MV/cm, are achieved.Since QDs' semiconductor lattice constant is larger than that in the barrier regions, a compressive stress in QDs and a tensile one inbarriers is induced in QDSL (Figure 2) with a consequently opposite sign for and along the z direction.Moreover, in the numerical evaluations of and , the spontaneous polarization provides a very low contribution with respect to the piezoelectricity.However, these values are too large to be considered negligible in the 1D Schrödinger equation along the z direction.Indeed, 1D minibands along the InGaN polar axis [0001] (Figure 3b) are split in a series of Wannier-Stark levels with envelope functions localized in a few dots, if the terms for and are included in the 1D Schrödinger equation.It is worth outlining that in the FEM simulations of Wannier-Stark levels, we have assumed fifteen periods in the z direction.This number of periods has been proven to be a good trade-off between the smaller computational domain and the solution accuracy.Indeed, our preventive investigations reveal that envelope functions tend to be localized in a number of 1D superlattice periods increasing with .In addition, numerical analysis obtained forcing the built-in electric field to zero reveals that in the samples analyzed, a few ( = 4) 1D minibands in the z direction are generally allowed.As a result, we can conclude that fifteen periods are enough to describe the spatial distribution of the wavefunction accurately.It is interesting to note that, sometimes, the number of QDs in which the envelope function is localized is the same as the order number , depending on the energy distance between consecutive Wannier-Stark levels given as: In particular, the abovementioned condition takes place for structures in which is larger than the QDs energy barrier as occurs for charge carriers in valence sub-bands in our samples.In Figure 4 energy levels for charge carriers in CB and in the A valence sub-band are shown.For both charge carriers, two types of envelope functions and then two groups of Wannier-Stark levels occur with an energy distance between consecutive allowed levels = 19.2meV.Although bands' edge However, these values are too large to be considered negligible in the 1D Schrödinger equation along the z direction.Indeed, 1D minibands along the InGaN polar axis [0001] (Figure 3b) are split in a series of Wannier-Stark levels with envelope functions localized in a few dots, if the terms for F z,dot and F z,bar are included in the 1D Schrödinger equation.
It is worth outlining that in the FEM simulations of Wannier-Stark levels, we have assumed fifteen periods in the z direction.This number of periods has been proven to be a good trade-off between the smaller computational domain and the solution accuracy.Indeed, our preventive investigations reveal that envelope functions ϕ t,epηq pzq tend to be localized in a number of 1D superlattice periods increasing witht.In addition, numerical analysis obtained forcing the built-in electric field to zero reveals that in the samples analyzed, a few (t max = 4) 1D minibands in the z direction are generally allowed.As a result, we can conclude that fifteen periods are enough to describe the spatial distribution of the wavefunction ϕ t,epηq pzq accurately.It is interesting to note that, sometimes, the number of QDs in which the envelope function ϕ t,epηq pzq is localized is the same as the order number t, depending on the energy distance between consecutive Wannier-Stark levels ∆E given as: In particular, the abovementioned condition takes place for structures in which ∆E is larger than the QDs energy barrier r V epηq z as occurs for charge carriers in valence sub-bands in our samples.
Photonics 2016, 3, 5 10 of 25 In Figure 4 energy levels for charge carriers in CB and in the A valence sub-band are shown.For both charge carriers, two types of envelope functions and then two groups of Wannier-Stark levels occur with an energy distance between consecutive allowed levels ∆E = 19.2meV.Although bands' edge bending in an ideal infinite superlattice leads to an infinite number of Wannier-Stark levels for each ϕ t,epηq pzq throughout the whole structure, in our 1D FEM simulation, we focus the investigation on Wannier-Stark levels related to localized envelope functions that are not negligible in one dot in the middle of the 1D FEM simulation domain.Then, five levels for each group are considered, and from the corresponding evaluation of the absorption coefficient, the information for the infinite structure can be extrapolated.In the sample of Figure 4, the groups of five energy levels for CB are separated, while for the A valence sub-band, they overlap.Often, nanostructures can be characterized by anisotropy, whichmodifies their responses to an external input as a function of its application direction.For the optical absorption, it is possible to estimate its dependence on light polarization.Our developed absorption coefficient model based on a Fermi's golden rule implementation done by using a numerical integration over the first Brillouin zone allows analysis of the light-matter interaction for different directions of q in order to have a more accurate investigation of light polarization influence.
Photonics 2016, 3, 5 10 of 25 bending in an ideal infinite superlattice leads to an infinite number of Wannier-Stark levels for each throughout the whole structure, in our 1D FEM simulation, we focus the investigation on Wannier-Stark levels related to localized envelope functions that are not negligible in one dot in the middle of the 1D FEM simulation domain.Then, five levels for each group are considered, and from the corresponding evaluation of the absorption coefficient, the information for the infinite structure can be extrapolated.In the sample of Figure 4, the groups of five energy levels for CB are separated, while for the A valence sub-band, they overlap.Often, nanostructures can be characterized by anisotropy,whichmodifies their responses to an external input as a function of its application direction.For the optical absorption, it is possible to estimate its dependence on light polarization.Our developed absorption coefficient model based on a Fermi's golden rule implementation done by using a numerical integration over the first Brillouin zone allows analysis of the light-matter interaction for different directions of in order to have a more accurate investigation of light polarization influence.Due to the different nature of A, B and C valence sub-band Bloch lattice functions, we have divided this analysis into two parts.At the first step, we have analyzed the influence on carriers' transition valence sub-band from A to CB minibands (Figure 5).Indeed, the A valence sub-band Bloch lattice functions, like heavy holes of ZB semiconductors, are characterized by only two directional and components, as is possible to see from Equation (4).On the contrary, B and C valence sub-bands' Bloch lattice functions present a directionality given by all three , and components.Furthermore, they are in a reciprocal relationship through the and parameters defined in Equation ( 5), depending on the In content of InGaN.
In  Due to the different nature of A, B and C valence sub-band Bloch lattice functions, we have divided this analysis into two parts.At the first step, we have analyzed the influence on carriers' transition valence sub-band from A to CB minibands (Figure 5).Indeed, the A valence sub-band Bloch lattice functions, like heavy holes of ZB semiconductors, are characterized by only two directional |Xy and |Yy components, as is possible to see from Equation (4).On the contrary, B and C valence sub-bands' Bloch lattice functions present a directionality given by all three |Xy, |Yy and |Zy components.Furthermore, they are in a reciprocal relationship through the a and b parameters defined in Equation ( 5), depending on the In content of InGaN.
In Figure 5, the absorption coefficients are shown as evaluated for the third sample in Table 2 for interminiband transitions from the A valence sub-band to CB minibands for (a) TE-and (b) TM-polarized incident light.The absorption coefficient profile in the two cases is similar, but values for TE-polarized light are almost three-times larger than for TM-polarized light.Generally, this is due to the InGaN anisotropy properties as a WZ semiconductor, but it is also influenced by strain, QD shape and In content in the structure.The highest peaks are obtained for transitions between the first minibands in the VB and CB.The maximum value of 273 cm ´1 occurs for interminiband transition from the 111 A valence sub-band to 112 CB for a photon energy of 3 eV (wavelength of 413 nm).The other higher peaks are related to interminiband transitions from 111 and 112 A valence sub-bands to 111 and 141 (with 411) CB minibands.Table 3 shows all interminiband absorption coefficient maximum peaks for the three QDSL structures described above and presented in Table 2.  3 shows all interminiband absorption coefficient maximum peaks for the three QDSL structures described above and presented in Table 2.  3).In Figure 6, the absorption coefficient spectra for interminiband transitions from the B and C valence sub-bands' to CB for TE-and TM-polarized light are shown.Furthermore, in these cases, absorption profiles are similar for different light polarizations, but absorption values for TM-polarized light are just over one and a half times greater than values for TE-polarized light for transitions involving minibands in the B valence sub-band, while for C, absorption values for TM-polarized light are a little less than one and a half times greater with respect to those for TE-polarized light.
The highest peaks are related to interminiband transitions from the B and C valence sub-band 121 (211) to 121 (211) CB minibands.For transitions from the B valence sub-band, an additional high peak is obtained for interminiband transitions from 121 (211) to 122 (212) in CB (Table 3).
Table 3. Interminiband absorption coefficient maximum peaks for In 0.4 Ga 0.6 N/GaN QDSL with dot sizes L β,dot = 4-6 nm, L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,bar = 2.5 nm evaluated with a Gaussian inhomogeneous linewidth broadening with a full width at half maximum (FWHM) broadening δ FW HM = 100 meV.The behavior of the interminiband absorption coefficient is presented in Figure 7 for all samples of Table 2 when QDs' base sizes are reduced, while keeping the other geometrical and compositional conditions the same.The behavior of the interminiband absorption coefficient is presented in Figure 7 for all samples of Table 2 when QDs' base sizes are reduced, while keeping the other geometrical and compositional conditions the same.The greater interaction of A valence sub-band minibands with TE-polarized light and of B and C with the TM light polarization described in previous paragraph is confirmed for all samples.Furthermore, for smaller QDs in asuperlattice, a minor number of minibands, corresponding to a thinner bandwidth of the absorption spectrum, is revealed.This is verified for all transitions involving minibands of different valence sub-bands, as is clearly visible in Figure 7 (a)-(c) (red and magenta contributions), and consequently, for the total absorption shown in Figure 7 (d).On the other hand, the situation is less regular for peaks values of the absorption coefficient.For transitions from the A valence sub-band minibands (a), smaller QDs lead to having higher peaks, while for transitions involving the B valence sub-band minibands (b), the trend is the opposite.In Figure 7 (c), a clear tendency is not visible, but this contribution to the total interminiband absorption coefficient The greater interaction of A valence sub-band minibands with TE-polarized light and of B and C with the TM light polarization described in previous paragraph is confirmed for all samples.Furthermore, for smaller QDs in asuperlattice, a minor number of minibands, corresponding to a thinner bandwidth of the absorption spectrum, is revealed.This is verified for all transitions involving minibands of different valence sub-bands, as is clearly visible in Figure 7a-c (red and magenta contributions), and consequently, for the total absorption shown in Figure 7d.On the other hand, the situation is less regular for peaks values of the absorption coefficient.For transitions from the A valence sub-band minibands (a), smaller QDs lead to having higher peaks, while for transitions involving the B valence sub-band minibands (b), the trend is the opposite.In Figure 7c, a clear tendency is not visible, but this contribution to the total interminiband absorption coefficient (d) is the most significant.However, in total, the interminiband absorption coefficient highest peaks of 2249 and 1613 cm ´1 for TM-and TE-polarized light, respectively, for a photon energy of 3.12 eV (wavelength of 397 nm) are obtained for the smallest QDs.These peaks are derived, above all, by the interminiband transition contributions from A, B and C valence sub-band to 111-112 CB minibands.

Interminiband Absorption Coefficient
Photonics 2016, 3, 5 13 of 25 (d) is the most significant.However, in total, the interminiband absorption coefficient highest peaks of 2249 and 1613 cm −1 for TM-and TE-polarized light, respectively, for a photon energy of 3.12 eV (wavelength of 397 nm) are obtained for the smallest QDs.These peaks are derived, above all, by the interminiband transition contributions from A, B and C valence sub-band to 111-112 CB minibands.The relatively small values of the absorption coefficient in Figures 3-5 is not surprising, since the exciton contribution has been forced to zero.Indeed, in these simulations, our intent is to investigate the influence of the light polarization and the different valence sub-bands on the absorption process.
All peak values and spectral widths of the absorption coefficient obtained through our developed model change by varying the spectral broadening for each transition introduced through a Gaussian inhomogeneous linewidth broadening ( ) (Equation (C2) in Appendix C).For the samples shown, a = 100 meV constant for each photon energy has been used.This model The relatively small values of the absorption coefficient in Figures 3-5 is not surprising, since the exciton contribution has been forced to zero.Indeed, in these simulations, our intent is to investigate the influence of the light polarization and the different valence sub-bands on the absorption process.
All peak values and spectral widths of the absorption coefficient obtained through our developed model change by varying the spectral broadening for each transition introduced through a Gaussian inhomogeneous linewidth broadening (δ FW HM ) (Equation (C2) in Appendix C).For the samples shown, a δ FW HM = 100 meV constant for each photon energy has been used.This model fitting parameter could be chosen with appropriate values variable with the photon energy in order to allow a better matching between the absorption coefficient predicted by the model and empirical measures.Indeed, linewidth broadening depends on geometrical (size of dots and regularity of their distribution) and compositional fluctuation (spatial variation of In content in the structure) that can influence in a different manner the absorption coefficient for various photon energies.The δ FW HM influence on the absorption coefficients is given in Figure 8, in which the exciton effect is also introduced for the second structure shown in Table 2.A broadening of Γ Broad = 7 meV for the exciton absorption is assumed as a good trade-off between very high exciton absorption coefficient peaks and its spectral width withrespect to that of the interminiband absorption obtained by the model.Furthermore, in this case, the value of broadening Γ Broad is a fitting parameter of the model for a more accurate matching with the actual empirical measures of the exciton absorption coefficient for the structure.
In Figure 8, the first highest three peaks visible on the blue curve for both TE-(a) and TM-polarized light (b) from left to right are related to the exciton effect for the B, A and C valence sub-band CB interminiband transitions, respectively.Even for exciton peaks, transitions involving minibands in the A valence sub-band are characterized by a greater interaction with TE light polarization, while for minibands in the C valence sub-band, the opposite situation occurs.As is visible in Figure 8, the exciton absorption coefficient peaks for transitions from the A and B valence sub-bands minibands are obtained for very near photon energy values (2.6 and 2.58 eV, corresponding to wavelengths of 477 and 480 nm, respectively, for the transition involving the first miniband in the A and B valence sub-bands).Then, both peak values are influenced by the other and, in this case, are higher for TE light polarization.When QDs' In concentration x changes in an In x Ga 1´x N/GaN QDSL, all of their electronic and optical properties are modified.For example, the barrier potentials for charge carriers in QDs can be increased or decreased according to the QD In x Ga 1´x N band gap variation with respect to the GaN of barrier regions.In Figure 9, interminiband absorption coefficients are shown for different QDSL samples with the same geometrical structure, but different QDs' In content x = 0.25, 0.325 and 0.4.Excitonic parameters are reported in Table 4.
Photonics 2016, 3, 5 14 of 25 fitting parameter could be chosen with appropriate values variable with the photon energy in order to allow a better matching between the absorption coefficient predicted by the model and empirical measures.Indeed, linewidth broadening depends on geometrical (size of dots and regularity of their distribution) and compositional fluctuation (spatial variation of In content in the structure) that can influence in a different manner the absorption coefficient for various photon energies.The influence on the absorption coefficients is given in Figure 8, in which the exciton effect is also introduced for the second structure shown in Table 2.A broadening of = 7 meV for the exciton absorption is assumed as a good trade-off between very high exciton absorption coefficient peaks and its spectral width withrespect to that of the interminiband absorption obtained by the model.Furthermore, in this case, the value of broadening is a fitting parameter of the model for a more accurate matching with the actual empirical measures of the exciton absorption coefficient for the structure.
In Figure 8, the first highest three peaks visible on the blue curve for both TE-(a) and TM-polarized light (b) from left to right are related to the exciton effect for the B, A and C valence sub-band CB interminiband transitions, respectively.Even for exciton peaks, transitions involving minibands in the A valence sub-band are characterized by a greater interaction with TE light polarization, while for minibands in the C valence sub-band, the opposite situation occurs.As is visible in Figure 8, the exciton absorption coefficient peaks for transitions from the A and B valence sub-bands minibands are obtained for very near photon energy values (2.6 and 2.58 eV, corresponding to wavelengths of 477 and 480 nm, respectively, for the transition involving the first miniband in the A and B valence sub-bands).Then, both peak values are influenced by the other and, in this case, are higher for TE light polarization.When QDs' In concentration x changes in an InxGa1-xN/GaN QDSL, all of their electronic and optical properties are modified.For example, the barrier potentials for charge carriers in QDs can be increased or decreased according to the QD InxGa1-xN band gap variation with respect to the GaN of barrier regions.In Figure 9, interminiband absorption coefficients are shown for different QDSL samples with the same geometrical structure, but different QDs' In content x = 0.25, 0.325 and 0.4.Excitonic parameters are reported in Table 4.In the structure with a greater x value, the band gap in QDs is lower, and the absorption spectrum is characterized by a red shift with respect to the samples with smaller x.For this reason, in In the structure with a greater x value, the band gap in QDs is lower, and the absorption spectrum is characterized by a red shift with respect to the samples with smaller x.For this reason, in Figure 9, the absorption spectrum of In 0.4 Ga 0.6 N/GaN QDSL (red line) extends starting from a photon energy of about 2.5 eV (a wavelength of about 500 nm) with respect to the structure with a lower x value (green and blue line) with a spectrum from about 2.7 and 2.8 eV.Further, the QDSL miniband number for x = 0.4 is greater than in other samples, so increasing its spectrum extension.In order to have absorption for greater wavelengths (lower photon energy), it is possible to increase the In concentration in the structure.Absorption coefficients for an In x Ga 1´x N/In 0.45 Ga 0.55 N QDSL for different x values are shown in Figure 10.Then, with x = 1 (InN/In 0.45 Ga 0.55 N QDSL),the absorption for a wavelength until about 1000 nm (around 1.2 eV) can be achieved with a spectrum until about 1.8 eV for interminiband transitions.This additional absorption contribution due to the miniband formation in superlattices increases the total absorption of an InN/In 0.45 Ga 0.55 N QDSL given by the other contribution related to interband transition into the barrier for high photon energy over 1.9 eV (about In 0.45 Ga 0.55 N band gap).Thus, this total absorption spectrum of an InN/In 0.45 Ga 0.55 N QDSL corresponds to that of an ideal intermediate band solar cell (IBSC) with one intermediate band and optimal band gaps values of 0.71, 1.24 and 1.95 eV for a maximum theoretical efficiency of this structure as presented in different works [53,54].Figure 9, the absorption spectrum of In0.4Ga0.6N/GaNQDSL (red line) extends starting from a photon energy of about 2.5 eV (a wavelength of about 500 nm) with respect to the structure with a lower x value (green and blue line) with a spectrum from about 2.7 and 2.8 eV.Further, the QDSL miniband number for x = 0.4 is greater than in other samples, so increasing its spectrum extension.In order to have absorption for greater wavelengths (lower photon energy), it is possible to increase the In concentration in the structure.Absorption coefficients for an InxGa1-xN/In0.45Ga0.55NQDSL for different x values are shown in Figure 10.Then, with x = 1 (InN/In0.45Ga0.55NQDSL),the absorption for a wavelength until about 1000 nm (around 1.2 eV) can be achieved with a spectrum until about 1.8 eV for interminiband transitions.This additional absorption contribution due to the miniband formation in superlattices increases the total absorption of an InN/In0.45Ga0.55NQDSL given by the other contribution related to interband transition into the barrier for high photon energy over 1.9 eV (about In0.45Ga0.55Nband gap).Thus, this total absorption spectrum of an InN/In0.45Ga0.55NQDSL corresponds to that of an ideal intermediate band solar cell (IBSC) with one intermediate band and optimal band gaps values of 0.71, 1.24 and 1.95 eV for a maximum theoretical efficiency of this structure as presented in different works [53,54].For lower x values (green and blue lines in Figure 10), charge carriers' barrier potential is lower and allows the formation of a smaller number of minibands with a spectrum around that of the exciton effect.Finally, Figures 6-8 reveal that the exciton effect increases the absorption coefficient of the order of 10 4 cm −1 depending on the light polarization and QD In content.Moreover, higher absorption values occur for TM light polarization.For lower x values (green and blue lines in Figure 10), charge carriers' barrier potential is lower and allows the formation of a smaller number of minibands with a spectrum around that of the exciton effect.Finally, Figures 6-8 reveal that the exciton effect increases the absorption coefficient of the order of 10 4 cm ´1 depending on the light polarization and QD In content.Moreover, higher absorption values occur for TM light polarization.
Thus, generally speaking, an increasing of the QD In content induces two main effects.The former is the energy bandgap reduction, resulting in the red shift absorption.The latter is the change of the piezoelectric effect as induced by an increased lattice mismatch between dot and barrier materials.Thus, generally speaking, an increasing of the QD In content induces two main effects.The former is the energy bandgap reduction, resulting in the red shift absorption.The latter is the change of the piezoelectric effect as induced by an increased lattice mismatch between dot and barrier materials.

Conclusions
In this paper, a model to evaluate miniband formation and the absorption coefficient for interminiband transitions of InGaN QDSL is presented.It introduces some simplifications as the cuboid shape for QDs and a wavefunction analysis based on spatial decoupling, but can consider the influence of strain, internal electric field due to spontaneous and piezoelectric polarizations of the material and can investigate the light polarization influence in an accurate manner.Furthermore,the excitonic effect related to the superlattice is taken into account.Through this model, the interminiband absorption coefficient of some samples has been evaluated showing very useful structures for many applications in which specific absorption spectra are needed (e.g., solar cells).Finally, some results show the InN/In0.45Ga0.55NQDSL to bea good candidate to make IBSCs.

A. Miniband Investigation
Investigation of QDSL minibands and charge carriers' steady states starts from one band single-particle model applied to electrons or holes in the A, B and C valence sub-bands.For each charge carrier in its own band, a 3D time-independent Schrödinger equation in the effective mass approximation has been considered as Equation (A1).We assume charge carriers of each band not interacting with those of other bands.Thus, the Schrödinger equation is: q q e pst q q e pst q q z e e where represents the reduced Plank constant and ê is the charge of an electron (with a negative sign) or hole (with positive sign).The functions indicate the potential energy.Moreover, is the effective mass tensor for electrons or holes in (A, B or C) valence sub-bands and evaluated by means of = , and being the electron (hole) relative effective mass (see Table 1) and free electron rest mass, respectively.

Conclusions
In this paper, a model to evaluate miniband formation and the absorption coefficient for interminiband transitions of InGaN QDSL is presented.It introduces some simplifications as the cuboid shape for QDs and a wavefunction analysis based on spatial decoupling, but can consider the influence of strain, internal electric field due to spontaneous and piezoelectric polarizations of the material and can investigate the light polarization influence in an accurate manner.Furthermore, the excitonic effect related to the superlattice is taken into account.Through this model, the interminiband absorption coefficient of some samples has been evaluated showing very useful structures for many applications in which specific absorption spectra are needed (e.g., solar cells).Finally, some results show the InN/In 0.45 Ga 0.55 N QDSL to bea good candidate to make IBSCs.

A. Miniband Investigation
Investigation of QDSL minibands and charge carriers' steady states starts from one band single-particle model applied to electrons or holes in the A, B and C valence sub-bands.For each charge carrier in its own band, a 3D time-independent Schrödinger equation in the effective mass approximation has been considered as Equation (A1).We assume charge carriers of each band not interacting with those of other bands.Thus, the Schrödinger equation is: ∇ 1 m ˚epηq px, y, zq `Vepηq px, y, zq ¯êF z pzq ff Ψ pst,epηq pq x ,q y q px, y, zq " E pst,epηq pq x ,q y q Ψ pst,epηq pq x ,q y q px, y, zq where represents the reduced Plank constant and ê is the charge of an electron (with a negative sign) or hole (with positive sign).The functions V epηq px, y, zq indicate the potential energy.Moreover, m ˚epηq px, y, zq is the effective mass tensor for electrons or holes in η (A, B or C) valence sub-bands and evaluated by means of m ˚epηq px, y, zq=m ˚epηq r px, y, zq m 0 , m ˚epηq r px, y, zq and m 0 being the electron (hole) relative effective mass (see Table 1) and free electron rest mass, respectively.Finally, the eigenfunctions Ψ pst,epηq pq x ,q y q px, y, zq and E pst,epηq pq x ,q y q represent the 3D envelope function and its eigenenergy for a given value of the vector q having components q x , q y and q z in the reciprocal space.Thus, in the envelope function paradigm, the steady-state wavefunction for electrons (holes) is given by the product of the 3D envelope function Ψ pst,epηq pq x ,q y q px, y, zq and the corresponding Bloch lattice function of semiconductor u CB 1{2 (u A 1{2 , u B 1{2 and u C 1{2 ).It is worth outlining that in Equation (A1), we have assumed the effective mass tensor as a function of the real-space coordinates, as a result of the influence of the heterostructure strain distribution.In addition, the spontaneous and piezoelectric polarizations are included in our modeling by means of the term F z pzq.Generally speaking, the built-in electric field must be considered in the vectorial form with all of its components spatially dependent.However, for the structures considered in this paper, the electric field component along the x and y directions can be considered negligible with respect to the zcomponent.Thus, in this sense, in Equation (A1), only the contributed F z pzq is considered.
In this context, it is worth noting than the inclusion of the strain, the spontaneous and piezoelectric polarizations represents a non-trivial generalization with respect to the models proposed in the literature [2][3][4][5][6][7], in which the electric field is neglected and the effective massesare assumed constant throughout QDs and barriers.In particular, we guess that the absence of the electric field along the z direction can represent a strong approximation, since spontaneous polarization, depending on temperature, and the piezoelectric effect, induced by technological features (see Section 2), are always present in InGaN devices, also in the absence of any heterostructure.At this step, if QDSL is characterized by a periodicity in the three orthogonal directions, Equation (A1) can be decoupled into three 1D equations describing 1D periodic structures [1][2][3][4][5][6][7].Then, a 3D envelope wavefunction Ψ pst,epηq pq x ,q y q px, y, zq is obtained as a product of three 1D envelope functions ϕ p,epηq q x pxq, ϕ s,epηq q y pyq and ϕ t,epηq pzq, as: Ψ pst,epηq pq x ,q y q px, y, zq " ϕ p,epηq q x pxq ¨ϕs,epηq q y pyq ¨ϕt,epηq pzq (A2) E pst,epηq pq x ,q y q " E p,epηq q x `Es,epηq where E p,epηq q x and E s,epηq q y represent the single particle carrier allowed energy levels of 1D minibands associated withthe periodic envelope functions ϕ p,epηq q x pxq and ϕ s,epηq q y pyq, like superlattice Bloch functions, as evaluated by using the Kronig-Penney model [55] applied along the x and y direction, respectively.On the contrary, E t,epηq and ϕ t,epηq pzq represent the Wannier-Stark energy levels and envelope functions along the z direction, respectively.Indeed, when the built-in electric field is applied in the z direction (WZ polar axis) of a superlattice structure, the translational invariance of its minibands' periodic envelope functions is broken.Minibands are split into groups of several single energy levels called a Wannier-Stark ladder, and their wavefunctions become localized in a few periods of the structure.In particular, if the 1D superlattice along the z direction suffers an electric field constant along the structure, the energy distance between the consecutive Wannier-Stark levels ∆E is equal to the product of the electric field and superlattice z-period.Moreover, the superscripts p and s represent the order numbers of 1D minibands along the x and y directions, respectively.Similarly, the t superscript is theorder number for the minibands along z when the built-in electric field is forced to zero.However, it is worth outlining that decoupling the 3D problem into three 1D ones can be mathematically possible if cuboid QDs' edges are orientated along the coordinates' directions (reference system) and if the electron or hole effective mass m ˚epηq px, y, zq tensor is assumed to be diagonal in the considered reference system.
In our studies, we consider samples with the QDs' growth direction corresponding to the InGaN polar crystallographic axis [0001] and the z direction of our reference system.Then, the electron or hole effective mass m ˚epηq px, y, zq tensor can be represented in the 1D approach by one of three scalar It is worth outlining that in this case, the computational domain must be chosen opportunely in order to describe the localization of the ϕ t,epηq pzq envelope functionsaccurately. 3D minibands are given by the composition of 1D minibands for two different directions obtained from Equation (A5) Wannier-Stark states in the z direction using Equation (A3) in each point (q x , q y , q z ) of the 3D superlattice reciprocal space.In this domain, for the axes' direction (generally labeled here with γ) corresponding to the related direction in real space, when 3D energy levels E pst,epηq γ occurred higher than r V epηq β , the continuum feature of InGaN in QDs and barriers isdominant, and optical properties can be assumed as those of the host (barrier regions) bulk material.

B. Piezo-Electro-Mechanical Modeling
As mentioned before, the generalization of our model is represented by the presence of the strain, spontaneous and piezoelectric polarizations.For this reason, it is worth briefly describingthe algorithm procedure used to integrate together the piezo-electro-mechanical modeling with the set of equations described in the previous section.In an In x Ga 1´x N/In y Ga 1´y N QDSL on a GaN buffer on a foreign substrate(see Section 2),there are two strains contributing.The first is related to lattice mismatch between the substrate with respect to the GaN buffer layer and the latter with respect to the In y Ga 1´y N host material.The second is related tothe lattice mismatch of QD heterostructures.However, the former can be neglected because its value could be reducible with an appropriate design of the structure or with specific technological innovation [26,39,40].The latter is intrinsic to the nature of QD heterostructures and is localized in a very small QD buried volume.Then, its effects have a great influence on the electric and optical properties of QDSL.Both miniband formation and optical absorption are modified through the electric field induced by the piezoelectric effect and by the direct influence of strain on the energy bands' edges and charge carrier effective masses.With the aim to realize self-consistent simulations, we have implemented an integrated algorithmic procedure based on home-made code and commercial software using a full-vectorial FEM.
In particular, for a given initial strain due to the lattice mismatch, the elastostatic linear continuum deformation of the heterostructure is evaluated as thefirst step by solving the following equations: where the subscripts i, j, l, m, kindicate the x, y or z directions for tensor components.In the case of the orthorhombic symmetry of cuboid QD, the initial normal strain satisfy the relationship ε 0 xx = ε 0 yy = ε a (biaxial strain) and ε 0 zz = ε c (principal diagonal elements of the initial strain tensor).In Equations (B1) and (B2), σ, ε and ε r represent the tensors of rank two of elastic stress, elastic strain and relative dielectric permittivity, respectively.Similarly, e F and C indicate the piezoelectric tensors of rank three and the elastic stiffness tensors of rank four, respectively.Moreover, D F and F represent the electric displacement field and electric field, respectively.The term P i is considered to be acting as a pyroelectric effect if a thermal change is applied or can represent the spontaneous polarization in wurtzite materials at a given temperature.In our context, P presents only the zcomponent, and the electric field F corresponds to that induced by the piezoelectric effect, without any external component.Finally, ε 0 is the vacuum dielectric permittivity.The energy levels in CB E CB σ and in all valence sub-bands E A σ , E B σ and E C σ with their relative effective masses m ˚CB r,zpβq , m ˚A r,zpβq , m ˚B r,zpβq and m ˚C r,zpβq can be found by applying the k¨p theory of Luttinger-Kohn and Bir-Pikus [8]: where the term V F represents the electric potential induced by electric field F.
It is worth outlining that the effective mass for electrons in CB is not explicitly evidenced since it is weakly influenced by the stress distribution.Moreover, in Equations (B3)-(B8), we have assumed ε xx = ε yy and ε xy = ε xz = ε yz = 0 as induced by the fact that the growth direction corresponds to the semiconductor polar axis, and the QDs are characterized by a cuboid shape with a square base.The coefficients a 1 and a 2 and D 1 , D 2 , D 3 , D 4 are WZ deformation potentials for CB and valence sub-bands, respectively, while A 1 , A 2 , A 3 , A 4 are WZ VB effective mass parameters.In InxGa1-xN materials, ε r , e F , C tensors components, WZ deformation potentials and VB effective masses depend on In content x as linear functions of those for InN and GaN values (Table 1).At this step, the spatial dependent variables , as well as the mean effective masses.Finally, the set of equations proposed in Appendix A is solved by the home-made code and 1D FEM Schrödinger module, allowing the minibands and the Wannier-Stark energy to be found.Then, the information about the electronic features of the QDSL is used to calculate the absorption spectrum as described in the following Appendix C.

C. Interminibands Absorption Coefficient
After evaluating the 3D minibands and charge carriers'steady-state envelope functions, the absorption coefficient for charge carriers' transitions from minibands into the A, B or C (η) valence sub-bands to those in CB in thermal equilibrium is calculated from a Fermi's golden rule implementation based on an integration into the first Brillouin zone (Equation (C1)), assuming a spectral Gaussian inhomogeneous linewidth broadening G (Equation (C2)) with an FWHM broadening δ FW HM depending onthe random structure sizes and compositional fluctuations as in [56].q ˇˇΨ r n η q Eˇˇˇ2 ˇˇ" r ¨Mη q ˇˇ2 G ´ ω ´Eg ´Er n e q ´Er n η q ¯dq x dq y dq z (C1) Variable energy value ω represents the photon energy; r is the light polarization direction; C 0 is defined in Equation (C3); and r N epηq and r n epηq are the total and progressive numbers, respectively, of allowed energy minibands for an electron (hole) in QDSL.Each integer progressive number r n epηq is associated withthe three integers pst labeling the 3D minibands (see Appendix A).In this work, the x (TE) and z (TM) directions are considered for light polarization r.M η q shown in Equation (C4) is the interband momentum matrix element for transition from valence sub-band η to CB, depending on the q vector of reciprocal space:

G
´ ω ´Eg ´Er n e q ´Er n η q ¯" 1 ?
2πσ c e ´p ω ´Eg ´Er n e q ´Er n η q q 2 2σ c 2 , σ c " δ FW HM 2 a 2lnp2q (C2) C 0 " π ê2 n r c 0 P 0 m 0 2 ω (C3) In Equation (C3), n r is the refractive index averaged on volume, and c 0 is the vacuum light velocity.The integral in Equation (C1) is made on the total volume of the first Brillouin zone Ω.For a number of structures, as in the bulk semiconductor, linewidth function G is replaced by a Dirac delta function.Moreover, the overlap integral ˇˇΨ r n e q ˇˇΨ r n η q ˇˇbetween charge carriers' envelope function and the interband momentum matrix element for a given light polarization ˇˇr ¨Mη q ˇˇcan be assumed q independent.Then, these become multiplicative factors of integration on the Dirac delta function from which, sometimes, the joint density of states is obtained using the parabolic model of q-space dispersions for the CB and VB [9].Instead, in our model for superlattices, in order to consider the non-parabolicity of minibands and the absorption coefficient light polarization dependence in a more accurate manner, each factor of the integrand is considered with its q dependence.It is numerically resolved as the sum of addends for different values of three q components q x , q y and q z related to a small q-space volume ∆Ω=∆q x ¨∆q y ¨∆q z around it and to its own energy value E g `Er n e q `Er n η q in which E r n e q and E r n η q are obtained by Equation (A3) in Appendix A. Actually, in our model, q z dependence is absent in Ψ r n epηq q and E r n epηq q due to the built-in electric field introduction along the z direction (see Appendix A).Only M η q is q z dependent, but ˇˇr ¨Mη q ˇˇ2 is averaged in this q-space direction.Inthis way, in the numerical solution of the integral in Equation (C1), a small q-space volume corresponds to ∆Ω =∆q x ¨∆q y ¨2π{T z in which T z is the 1D superlattice spatial period along the z direction.The M η q evaluation needs the superlattice Kane parameters.They are estimated as volumetric mean values between those in QDs and barrier regions of their values K η β in Equation (C5) and K η z in Equation (C6) [10,34] defined in the k¨p formalism.Generally, Kane parameters are assumed strain independent.Instead, as donein [11], we have considered the crystal deformation influence.Then, we introduce in our model a mean charge relative reduced mass value m ˚eη r,βpzq strain dependent on Kane parameters formulas (all other semiconductor and particles variables are defined previously):

D. Exciton Absorption
In MQDs' structure in which each nanostructure can be assumed as a single QD (the absence of wavefunction overlap among those of different QDs), the main absorption coefficient component is related to QDs' bound excitons [35,43,47,[49][50][51][52].Instead, in this paper, our investigation is focused on absorption coefficient components derived by thesuperlattice nature of high density QD heterostructures.For this reason, in this Appendix, we estimate the contribution due to the free exciton derived from the charge Coulomb interaction between a single electron in the first minibands of the CB and valence sub-bands.For each of these excitons, a Bohr radius a Bohr , a Rydberg constant R ẙ and a reduced electron-hole effective mass m The contribution of the absorption coefficient of this exciton effect (α η exc ) based on the Elliott formula and assuming a spectral Lorentzian homogeneous linewidth broadening is shown in Equation (D3) [57,58]:  For the investigation of free excitons in Equation (D3), envelope functions Ψ e and Ψ η for the first minibands of the CB and η valence sub-band, respectively, M η and charge carrier energy levels E e and E η are considered for q = 0.For this reason, M η has been averaged over the solid angle in real space as in the typical analysis for an isotropic bulk semiconductorand, so, assumes a value of one-third of the Kane parameters.Anisotropic QDSLs' response to different light polarization occurs assuming the different values of K η β and K η z shown in Equations (C5) and (C6) for TE-and TM-polarized light, respectively, as done for interminiband absorption (see Appendix C).

Figure 1 .
Figure 1.In x Ga 1´x N/In y Ga 1´y N quantum dot superlattice (QDSL) 2D schematic view.Each cuboid In x Ga 1´x N QD (green) with sizes L x,dot = L y,dot = L β,dot (cuboid QD with square base) and L z,dot is embedded in an In y Ga 1´y N matrix (x > y) with inter-dot distances L x,bar = L y,bar =L β,bar and L z,bar .

Figure 3 .
Figure 3.In 0.4 Ga 0.6 N/GaN QDSL electron minibands for dot sizes L β,dot = 6 nm and L z,dot = 3 nm and inter-dot distances of L β,bar = 5 nm and L z,bar = 2.5 nm along the (a) [[100]] and (b) [[001]] directions in reciprocal q-space.The red line represents the electron energy potential in the barriers' region.

Figure 4 .
Figure 4. Wannier-Stark ladder energy levels (red lines for charge carriers in CB and black lines in the A valence sub-band) along the polar axis direction in an In0.4Ga0.6N/GaNQDSL with dot sizes = 6 nm, = 3 nm and inter-dot distances = 5 nm and = 2.5 nm.The subplot insets are the enlargement for (a) CB and (b) the A valence sub-band.Arbitrarily, an energy value of 0 eV has been assigned to the minimum of the CB edge into the dot in the first period shown.

Figure 5 ,
the absorption coefficients are shown as evaluated for the third sample in Table2for interminiband transitions from the A valence sub-band to CB minibands for (a) TE-and (b) TM-polarized incident light.The absorption coefficient profile in the two cases is similar, but values for TE-polarized light are almost three-times larger than for TM-polarized light.Generally, this is due to the InGaN anisotropy properties as a WZ semiconductor, but it is also influenced by strain, QD shape and In content in the structure.The highest peaks are obtained for transitions between the first minibands in the VB and CB.The maximum value of 273 cm −1 occurs for

Figure 4 .
Figure 4. Wannier-Stark ladder energy levels (red lines for charge carriers in CB and black lines in the A valence sub-band) along the polar axis direction in an In 0.4 Ga 0.6 N/GaN QDSL with dot sizes L β,dot = 6 nm, L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,bar = 2.5 nm.The subplot insets are the enlargement for (a) CB and (b) the A valence sub-band.Arbitrarily, an energy value of 0 eV has been assigned to the minimum of the CB edge into the dot in the first period shown.
Photonics 2016, 3, 5 11 of 25 interminiband transition from the 111 A valence sub-band to 112 CB for a photon energy of 3 eV (wavelength of 413 nm).The other higher peaks are related to interminiband transitions from 111 and 112 A valence sub-bands to 111 and 141 (with 411) CB minibands.Table3shows all interminiband absorption coefficient maximum peaks for the three QDSL structures described above and presented in Table2.

Figure 5 .
Figure 5. Absorption coefficient for interminiband transition from A valence sub-band to CB minibands for (a) TE-and (b) TM-polarized light in In0.4Ga0.6N/GaNQDSL with dot sizes = 6 nm, = 3 nm and inter-dot distances = 5 nm and = 2.5 nm.Integer numbers near the highest peaks specify the CB miniband involved in the absorption process.In Figure6, the absorption coefficient spectra for interminiband transitions from the B and C valence sub-bands' to CB for TE-and TM-polarized light are shown.Furthermore, in these cases, absorption profiles are similar for different light polarizations, but absorption values for TM-polarized light are just over one and a half times greater than values for TE-polarized light for transitions involving minibands in the B valence sub-band, while for C, absorption values for TM-polarized light are a little less than one and a half times greater with respect to those for TE-polarized light.The highest peaks are related to interminiband transitions from the B and C valence sub-band 121 (211) to 121 (211) CB minibands.For transitions from the B valence sub-band, an additional high peak is obtained for interminiband transitions from 121 (211) to 122 (212) in CB (Table3).

Figure 5 .
Figure 5. Absorption coefficient for interminiband transition from A valence sub-band to CB minibands for (a) TE-and (b) TM-polarized light in In 0.4 Ga 0.6 N/GaN QDSL with dot sizes L β,dot = 6 nm, L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,bar = 2.5 nm.Integer numbers near the highest peaks specify the CB miniband involved in the absorption process.

Figure 6 .
Figure 6.Absorption coefficient for interminiband transitions from the B valence sub-band to CB for (a) TE-and (b) TM-polarized light and interminiband transitions from C valence sub-band to CB minibands for (c) TE-and (d) TM-polarized light in In0.4Ga0.6N/GaNQDSL with dot sizes = 6 nm, = 3 nm and inter-dot distances = 5 nm and = 2.5 nm.Integer numbers near the highest peaks specify the CB minibandinvolvedin the absorption process.

Figure 6 .
Figure 6.Absorption coefficient for interminiband transitions from the B valence sub-band to CB for (a) TE-and (b) TM-polarized light and interminiband transitions from C valence sub-band to CB minibands for (c) TE-and (d) TM-polarized light in In 0.4 Ga 0.6 N/GaN QDSL with dot sizes L β,dot = 6 nm, L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,dot = 2.5 nm.Integer numbers near the highest peaks specify the CB minibandinvolvedin the absorption process.

Figure 7 .
Figure 7. Absorption coefficient for the (a) A, (b) B, (c) C and (d) total valence sub-band CB interminiband transitions for TE-and TM-polarized light in In 0.4 Ga 0.6 N/GaN QDSL with variable dot sizes L β,dot = 4-6 nm, L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,bar = 2.5 nm.

Figure 8 .
Figure 8.Total interminiband transition absorption coefficient with theexcitonic effect for (a) TE-and (b) TM-polarized light in In 0.4 Ga 0.6 N/GaN QDSL with dot sizes L β,dot = 5 nm, L z,dot = 3 nm and inter-dot distances L β,bar = 5 nm and L z,bar = 2.5 nm evaluated with a spectral broadening of Γ Broad = 7 meV and δ FW HM = 35-100 meV for the excitonic and interminibands' absorption contributions, respectively.

Figure 9 .
Figure 9.Total interminiband transition absorption coefficient with superlattice excitonic effect for (a) TE-and (b) TM-polarized light in In x Ga 1´x N/GaN QDSL with dot sizes L β,dot = 5 nm, L z,dot = 3 nm inter-dot distances L β,bar = 5 nm, L z,bar = 2.5 nm and QD In content x = 0.25-0.4.

Figure 10 .
Figure 10.Total interminiband transition absorption coefficient with the superlattice excitonic effect for (a) TE-and (b) TM-polarized light in In x Ga 1´x N/In 0.45 Ga 0.55 N QDSL with dot sizes L β,dot = 5 nm, L z,dot = 3 nm, inter-dot distances L β,bar = 5 nm, L z,bar = 2.5 nm and QDs In content x = 0.65-1.

2 ;
β labels one of the x or y directions and m ˚epηq z pzq.Each of these effective masses assumes different values in QDs with respect to the barrier regions.Thus, in the context of the decoupling approach, we assume the 1D periodic electron (hole) energy barriers calculated by means of an opportune averaging procedure in order to take into account the strain, spontaneous and piezoelectric polarizations.According to theKronig-Penney model in the x and y directions, the following dispersion equationfor pβ " x, yq(A6) where T β " L β,bar `Lβ,dot indicates the 1D superlattice period.The terms m epηq β,bar and m epηq β,dot represent the average values of m ˚epηq β pβq for the electron (hole) into barriers and QDs, respectively.Allowed energy levels for 1D minibands alongthe x and y directions are obtained as a function of q β varying values of E epηq β in the dispersion Equation (A5).To complete the 1D modeling, the Wannier-Stark energy levels E t,epηq and the envelope functions ϕ t,epηq pzq in the z direction are evaluated through FEM simulations of the 1D Schrödinger equation for fixed values of r and F z pzq.
zpβq , m ˚A r,zpβq , m ˚B r,zpβq and m ˚C r,zpβq are averaged into the volume of central dots (dot values) and of the barrier around it (barrier values) into the computational domain for the x, y and z directions, in order to evaluate the potential energy barrier's mean values r V

m ˚e x r m ˚e y r m ˚e z `1 3
˚eη exc can be evaluated by means of Equation (D1):a Bohr " 4πε r ε b r m ˚η x r m ˚η y r m ˚η z ‚´1 (D1)The values r m ˚e ν and r m ˚η ν are the semiconductor effective masses for electrons and holes in the x, y and z directions (ν), averaged on a period T ν :