Analytical Expressions for Numerical Characterization of Semiconductors per Comparison with Luminescence

Luminescence is one of the most important characterisation tools of semiconductor materials and devices. Recently, a very efficient analytical set of equations has been applied to explain optical properties of dilute semiconductor materials, with an emphasis on the evolution of peak luminescence gain with temperature and its relation to sample quality. This paper summarizes important steps of the derivation of these expressions that have not been presented before and delivers a theoretical framework that can used to apply exactly solvable Hamiltonians for realistic studies of luminescence in various systems.


Introduction
Materials development requires characterization techniques. Among them, photoluminescence, or in more general terms, radiation emission due to different excitation mechanisms [1], is a very powerful tool to study semiconductor materials and map specific characteristics of new devices [2,3] for applications from the THz-Mid Infrared (TERA-MIR) to ultraviolet ranges [4][5][6]. A one-to-one correspondence between measured spectra and fundamental materials properties requires a clear theoretical model, ideally easy to understand and to programme, but at the same time with microscopic information for conclusive interpretation and as free as possible from phenomenological parameters.
A recent theoretical effort led to the development of analytical solutions for the interband polarization, which plays the selfenergy role in the Dyson equation for the Photon Green's functions [7], which have been applied them to study photoluminescence of Coulomb-correlated semiconductor materials. The accuracy of the resulting easily programmable solutions has been demonstrated by consistently explaining the low temperature s-shape of the luminescence peak of dilute semiconductors, such as ternary GaAsSb, InAsN, and quaternary InAs(N,Sb) [7][8][9]. The interplay of homogeneous versus inhomogeneous broadening at low and high temperatures are described, together with the relevance of many body effects, which are in very good agreement with experiments [10][11][12][13]. A similar set of equations was also used to study nonlinearities in GaAs-AlGasAs and GaAsN-AlGaAs superlattices [14][15][16]. The superlattice case is particularly noteworthy for room temperature GHz nonlinear multiplication into the THz range [17,18]. This paper has two objectives: to show hitherto unpublished details of the mathematical steps that lead to the equations used in Refs. [7][8][9][14][15][16], and to draw a bridge between the luminescence and nonlinear absorption calculations in superlattices. It is organized as follows: the main steps involving manipulations of hypergeometric functions that characterize the solution of the Hulthén potential problem are delivered. Next, a direct connection between luminescence and absorption equations is given together with a connection with generalized semiconductor Bloch equations, with potential for the study of polaritons in superlattices within a dielectric approach such as that used to investigate valence band THz polaritons and antipolaritons [6]. A brief conclusion follows.

Mathematical Model
In order to make this paper self-contained, some of the steps shown in Ref. [7] are followed to guide the presentation towards the more complete derivation presented here.

Integro-Differential Equation for the Power Spectrum
Luminescence, or equivalently the optical power density spectrum I(ω), is described quantum mechanically by the Poynting vector, which can be directly related to the transverse polarization function P, which is the selfenergy in the Dyson equations for the transverse photon Green's function components in the Keldysh formalism [7,19,20]. All the quantities presented in this paper are considered in frequency space, i.e., evaluated at steady state.
The free photon Green's function represents the photons propagating without any interaction with the medium. When carriers are injected the transverse polarization function P, which is the selfenergy in photon Green's function Dyson equation, determines how the excited medium modifies the photon propagation. The lesser Keldysh component P < is proportional to the carriers recombination rate and yields the number of emitted photons per unit area. It thus governs the power emission spectrum, as seen in Equation (1). The imaginary and real parts of P r are, respectively, proportional to absorption and gain and refractive index changes, since the dielectric function of the medium reads (ω) = 1 − c 2 ω 2 P r (ω) as shown in Ref. [19]. The starting points for the results derived here are the equations in Refs. [7,19,20].
Here, e, c, Ω, and Π denote, respectively, the electron charge, the speed of light, the sample volume and the velocity matrix element, which is the expectation value of the velocity operator, i.e., the momentum operator divided by the electron mass. It stems directly from the fact that current is charge times velocity. The formal definition of the transverse polarization function selfenergy in terms of functional derivatives is P = − 4π c δJ δA , where J and A are, respectively, expectation values of the induced current and vector potential operators. The full expression involves labels along the Keldysh contour and is tensorial. A complete discussion is beyond the scope of this paper. For details see Refs. [19,20].
The crystal momentum → k is a consequence of Fourier transforming from real space. Likewise, ω and thus the photon energy ω stem from a corresponding Fourier transformation from time to the frequency domain. The matrix element satisfies the integro-differential equation [20], where W is the screened Hulthén potential [21][22][23][24]. Furthermore, Electrons or holes are labelled, respectively, by λ = {e, h}, the renormalized energies e λ ( → k ), and dephasing Γ λ are calculated from the real and imaginary parts of the selfenergy in the Dyson equation for the retarded carriers Green's functions. This paper focuses on quasi-equilibrium luminescence and on three dimensional (bulk semiconductors) with one conduction and one valence band.
Under these conditions, f λ denotes a Fermi function characterized by a chemical potential µ λ and the spectral function in Equation (4) for each particle, derived from components of the carriers' Green's function in the Keldysh formalism, readŝ The next step is to re-write the last term in Equation (4) by means of the identity (6) and to approximate this factor by where µ = µ e + µ h is the total chemical potential, where T is the temperature in Kelvins and K B is the Boltzmann constant. Different versions of this approximation has been used before in phenomenological approaches for absorption Refs. [21][22][23] and delivered good agreement with experiments (see details and further references in Ref. [23]). Within the Keldysh Green's functions, context, a detailed derivation of its application is given in Ref. [20]. The fully numerical solutions of the equations that use this version of the approximation have given very good agreement with both single beam and pump-probe luminescence [20,25]. Its usefulness has been further confirmed recently by the good agreement between the analytical solutions shown here and the experimental luminescence of dilute semiconductors [7][8][9]14].
is an excellent approximation. Note that this theory is applied for photon energies around the semiconductor bandgap and thus in Figures 1 and 2, Under these conditions, denotes a Fermi function characterized by a chemical potential and the spectral function in Equation (4) for each particle, derived from components of the carriers' Green's function in the Keldysh formalism, reads The next step is to re-write the last term in Equation (4) by means of the identity and to approximate this factor by 1 − ( ) − ( − ) ≈ tanh (ℏ ) , where = + is the total chemical potential, where T is the temperature in Kelvins and is the Boltzmann constant. Different versions of this approximation has been used before in phenomenological approaches for absorption Refs. [21][22][23] and delivered good agreement with experiments (see details and further references in Ref. [23]). Within the Keldysh Green's functions, context, a detailed derivation of its application is given in Ref. [20]. The fully numerical solutions of the equations that use this version of the approximation have given very good agreement with both single beam and pump-probe luminescence [20,25]. Its usefulness has been further confirmed recently by the good agreement between the analytical solutions shown here and the experimental luminescence of dilute semiconductors [7][8][9]14].   ( , ), the term in curly braces in Equation (6) which has been approximated by ( , ) ≈ 1, evaluated at ℏ = for bulk GaAs at low temperatures, where the dephasing is typically small Γ . In this case, only a small range of detunings − ′ contribute to the integral in Equation (4). Thus the range chosen in the x-axis, from zero to approximately twice the exciton binding energy (2 ), is even larger than necessary. (a) T = 10 K; (b) T = 20 K.
Low temperature luminescence is typically performed with a small density of injected carriers. Very good agreement of this theory with results from different experimental teams for a variety of materials has been obtained with carrier densities around 10 15 carriers/cm 3 [7][8][9], further justifying the range of densities in the y-axis. The theory has also been used for high temperatures and high densities to investigate optical nonlinearites [14][15][16], and this range is illustrated in Figure 2.  (6) which has been approximated by FF(ω, ω ) ≈ 1, evaluated at ω = E g for bulk GaAs at low temperatures, where the dephasing is typically small Γ λ . In this case, only a small range of detunings ω − ω contribute to the integral in Equation (4). Thus the range chosen in the x-axis, from zero to approximately twice the exciton binding energy (2e 0 ), is even larger than necessary. (a) T = 10 K; (b) T = 20 K.
Low temperature luminescence is typically performed with a small density of injected carriers. Very good agreement of this theory with results from different experimental teams for a variety of materials has been obtained with carrier densities around 10 15 carriers/cm 3 [7][8][9], further justifying the range of densities in the y-axis. The theory has also been used for high temperatures and high densities to investigate optical nonlinearites [14][15][16], and this range is illustrated in Figure 2. ( , ), the term in curly braces in Equation (6) which has been approximated by ( , ) ≈ 1, evaluated at ℏ = for bulk GaAs at higher temperatures and high densities, where the dephasing Γ is larger than in the low temperature case. In this case a wider range of detunings − ′ contribute to the integral in Equation (4). Thus, the range chosen in the x-axis is even larger Relevant dephasing mechanisms such as electron-electron, electron-phonon and electronimpurity scattering can be added to the selfenergy [17,18], and the resulting Γ is frequency and momentum dependent. However, in what follows, it is replaced by averaged values, leading to a simple approximation for { ( , } consistent with the Ansatz solution, where Γ = Γ + Γ and ≡ tanh (ℏ − )/2 . In 3D, the material resonance energy is: Δ = ℏ + , where 1 * = 1 + 1 . The bandgap is given by the sum of the fundamental band gap , and a many body renormalisation term Δ where , denote, respectively the electron and hole effective masses. The equation for , simplifies to: The total dephasing Γ will determine the luminescence linewidth. Thus, it can be treated as a phenomenological parameter used to interpret data, and at the same time estimate the strength of the scattering and dephasing processes [7][8][9] by comparison of adjusted data with microscopic calculations derived from the relevant selfenergies [17,18]. At this point, the Kubo-Martin-Schwinger (KMS) relation under the form derived in Ref. [20] can be applied to Equation (8), together with the auxiliary variable: Λ , = , (ℏ )/( ) , leading to the relation: Expressing , from Equation (8) in terms of Λ , the corresponding integro-differential equation becomes  (6) which has been approximated by FF(ω, ω ) ≈ 1, evaluated at ω = E g for bulk GaAs at higher temperatures and high densities, where the dephasing Γ λ is larger than in the low temperature case. In this case a wider range of detunings ω − ω contribute to the integral in Equation (4). Thus, the range chosen in the x-axis is even larger than necessary. (a) T = 150 K, (b) T = 300 K.
Relevant dephasing mechanisms such as electron-electron, electron-phonon and electron-impurity scattering can be added to the selfenergy [17,18], and the resulting Γ λ is frequency and momentum dependent. However, in what follows, it is replaced by averaged values, leading to a simple approximation for Im{P 0 r (k, ω)} consistent with the Ansatz solution, where Γ = Γ e + Γ h and ϑ ≡ tanh[β( ω − µ)/2]. In 3D, the material resonance energy is: The bandgap E g is given by the sum of the fundamental band gap E 0 g , and a many body renormalisation term ∆E g where m e , m h denote, respectively the electron and hole effective masses. The equation for P r ( → k , ω) simplifies to: The total dephasing Γ will determine the luminescence linewidth. Thus, it can be treated as a phenomenological parameter used to interpret data, and at the same time estimate the strength of the scattering and dephasing processes [7][8][9] by comparison of adjusted data with microscopic calculations derived from the relevant selfenergies [17,18]. At this point, the Kubo-Martin-Schwinger (KMS) relation under the form derived in Ref. [20] can be applied to Equation (8), together with the auxiliary variable: Λ( , leading to the relation: . Before proceeding, the Hulthén potential [21][22][23][24] should be revised. The usual approximation for a static 3D screened potential is the Yukawa potential, W Y ( , has known analytical solutions that have proven to be very useful for the description of bulk absorption [22]. Recent applications have confirmed its relevance to explain experimental luminescence studies [7][8][9]. Figure 3 shows that, in the range of carrier densities and temperatures of interest, the Yukawa potential can be replaced by the Hulthén potential with negligible differences in numerical values. ) . Before proceeding, the Hulthén potential [21][22][23][24] should be revised. The usual approximation for a static 3D screened potential is the Yukawa potential, However, the corresponding Schrödinger equation does not have known analytical solutions. In contrast, the Hulthén potential: (| |) = −2 /((exp (2 | |) − 1)), has known analytical solutions that have proven to be very useful for the description of bulk absorption [22]. Recent applications have confirmed its relevance to explain experimental luminescence studies [7][8][9]. Figure 3 shows that, in the range of carrier densities and temperatures of interest, the Yukawa potential can be replaced by the Hulthén potential with negligible differences in numerical values. Both cases depend on the temperature T and carrier density N through the inverse screening length and the inset explains the results, because increases with increasing carrier density and with decreasing temperature. The dot-dashed (cyan) curve is for T = 10 K, while the double-dot-dashed (orange) curve is for T = 300 K.
The Fourier transform of the Hulthén potential has an analytical expression, where Ω is the sample volume, is the Trigamma function [26], is the screening wavenumber and by including in at Equation (11), = / . Analytical approximations for and are given in Ref. [7]. Note that the bandgap renormalization including Coulomb hole and screened exchange corrections reads Both cases depend on the temperature T and carrier density N through the inverse screening length κ and the inset explains the results, because κ increases with increasing carrier density and with decreasing temperature. The dot-dashed (cyan) curve is for T = 10 K, while the double-dot-dashed (orange) curve is for T = 300 K.
The Fourier transform of the Hulthén potential has an analytical expression, where Ω is the sample volume, ψ is the Trigamma function [26], κ is the screening wavenumber and by including ϑ in W at Equation (11), 0 = 0 /ϑ. Analytical approximations for µ and κ are given in Ref. [7]. Note that the bandgap renormalization including Coulomb hole and screened exchange corrections reads The Fermi functions f e , f h are evaluated at the peak of the spectral function for each particle, i.e., in Equation (5), ω = e λ ( → k ). More details are given in Ref. [7]. Equation (13) goes beyond phenomenological term for the bandgap shift [21][22][23], and also, in contrast to those, here we can in principle take into account a reduction in the Coulomb interaction due to phase space filling through the factor ϑ. Note however that in the range of carrier densities and temperatures of interest ϑ ≈ 1 i.e., 0 ≈ 0 , as shown in Figure 4.  At quasi-equilibrium, used in Refs. [7][8][9][14][15][16]22], the total chemical potential is calculated self-consistently with the many body renormalization of the bandgap and can be written exactly as = + , where is the total free carrier chemical potential calculated from the bottom of each band. In other words, the inversion factor can be equivalently written as ≡ tanh (ℏ − − )/2 , and it is now clear why Figure 4 has the detuning ℏ − in the x-axis. Furthermore, the 3D exciton binding energy for GaAs is 4.2 meV and there is no luminescence of absorption for a detuning below 4.2 meV, unless there are deep levels due to impurities and defects, which are not considered here. Thus the approximation, ≈ 1 i.e., ≈ for the dielectric constant used in the Hulthén potential is clearly excellent in the low power luminescence case. Nonlinear absorption studies are only meaningful away from population inversion leading to optical gain, i.e., the studies are in the range ≥ 1. Thus a decreasing occupation reflects phase space feeling and even for ≠ the approach is valid. In order to study the gain regime, the approximation used in the literature is to make at the Hulthén potential = and consider the inversion factor only on the right hand side of Equation (8). In other words, in the traditional "plasma theories" for bulk semiconductor absorption and gain, phase space filling ( ≠ 1) is not taken into account. Since in the high density case where gain develops, the Hulthén potential decreases due to screening, which described by large in Equation (12), there is still good agreement with experiments. See e.g., Refs. [22,23]. Equation (11) can now be Fourier-transformed Here Ω is the sample volume, ( ) denotes the Dirac delta function. Expanding Λ( , ) in the basis of eigenstates of the Hamiltonian: ℋ = − ℏ ∇ − ( ), At quasi-equilibrium, used in Refs. [7][8][9][14][15][16]22], the total chemical potential µ is calculated self-consistently with the many body renormalization of the bandgap E g and can be written exactly as where µ is the total free carrier chemical potential calculated from the bottom of each band. In other words, the inversion factor can be equivalently written as ϑ ≡ tanh β( ω − E g − µ)/2 , and it is now clear why Figure 4 has the detuning ω − E g in the x-axis. Furthermore, the 3D exciton binding energy for GaAs is 4.2 meV and there is no luminescence of absorption for a detuning below 4.2 meV, unless there are deep levels due to impurities and defects, which are not considered here. Thus the approximation, ϑ ≈ 1 i.e., 0 ≈ 0 for the dielectric constant used in the Hulthén potential is clearly excellent in the low power luminescence case. Nonlinear absorption studies are only meaningful away from population inversion leading to optical gain, i.e., the studies are in the range ϑ ≥ 1. Thus a decreasing occupation reflects phase space feeling and even for 0 = 0 the approach is valid. In order to study the gain regime, the approximation used in the literature is to make at the Hulthén potential 0 = 0 and consider the inversion factor only on the right hand side of Equation (8). In other words, in the traditional "plasma theories" for bulk semiconductor absorption and gain, phase space filling (ϑ = 1) is not taken into account. Since in the high density case where gain develops, the Hulthén potential decreases due to screening, which described by large κ in Equation (12), there is still good agreement with experiments. See e.g., Refs. [22,23].
Equation (11) can now be Fourier-transformed Here Ω is the sample volume, δ( → r ) denotes the Dirac delta function. Expanding Λ( → r , ω) in the basis of eigenstates of the Hamiltonian: Thus, Equation (15) can be rewritten as Projection onto state ν yields Substitution into Equation (16) Λ Fourier-transforming back to k-space From , a closed expression can be obtained.
where a factor 2 for spin has been explicitly written out of the summation over all quantum numbers. Introducing and combining Equations (1), (2), (10) and (21) leads to where δ Γ = 1 π Γ ( ω−E g −E ν ) 2 +Γ 2 reduces to a Dirac delta function for Γ → 0. The velocity matrix element is expressed in terms of the dipole moment matrix element and the fundamental bandgap as |Π| = (E 0 g / )|S|x|X|. Next, the Schrödinger Equation for the Hulthén potential must be solved, so that ψ ν (0) can be inserted in Equation (23). The first step is to separate the wavefunction in radial and angular parts. The label ν thus spans the set {n, l, m}, The corresponding Schrödinger Equation, which is a generalized Wannier equation [23] can be cast in the form: The energy eigenvalues depend only on the {n, l} quantum numbers, and thus we can replace E n by E nl . Introducing the 3D Rydberg e 0 and Bohr Radius a 0 , as well λ = 2κ, g = 1/(κa 0 ) and Note that the angular momentum operator has been applied to the wavefunction directly from Equation (25) to Equation (26), i.e., L 2 ψ ν ( → r ) = l(l + 1) 2 ψ ν ( → r ). Only solutions that do not vanish at r = 0 contribute to the emitted power, so l = 0 is selected. The labels "nl" will be dropped at the moment to simplify the notation, Introducing u = r f , and β = − λ 2 , Equation (27) is transformed into The auxiliary variables, z = 1 − e −λr and w = u z(1−z) β , lead to the equation which reduces to the Hypergeometric Equation [26,27], The generalized Wannier Equation, Equation (25), has two types of solutions: bound states for ν < 0 and unbound solutions for ν > 0. The wavefunctions and eigenvalues are thus different and it makes sense to study each case separately and then add all contributions when a sum over all possible ν as required from Equation (23).

Continuum States
The unbound solutions that make a continuum have positive eigenvalues, ν > 0, and thus imaginary β ν = i √ ν /λ. Dropping labels to simplify the development in the next equations yields which can be written for simplicity as Note that the transformation z = 1 − e −λr is being used. The solution that will be later inserted in Equation (23), will be normalized in a sphere of radius r = R and asymptotic solutions, obtained a large radius R will be investigated. Next, Equation (15.3.6) from Ref. [26] is used, i.e., F(a, b, c; z) = . Furthermore, note that lim R→∞ z = 1 and for all values of ξ, λ, δ, F(ξ, λ, δ; 0) = 1. Thus,  [26]), gives in the asymptotic limit Leading to asymptotic forms of u, |u| 2 The normalization constant is thus given by (52) However, |Γ(−2i|β|| 2 = π 2|β|sinh(2π|β|) , see e.g., Equations (6.1.29) and (6.1.31) of Ref. [26], plus a little algebra deliver the continuum normalization constant The required value of the wave function at the origin can thus be expressed as Next, note that lim r→0 z = 0, but lim r→0 z r = λ, leading to The sum of continuum states becomes an integral, ∑ ν . . . = Rλ π ∞ 0 . . . Introducing 2 for spin and changing variables, the continuum contribution becomes which combines with the bound states to deliver the power spectrum Here, ζ = ω−E g e 0 , e n = − 1 n 2 (1 − n 2 g ) 2 , I 0 = e 2 ω 2 |Π| 2 πe 0 c 3 a 3 0 and the square of the velocity matrix element where the spin orbit shift, the free-carrier bandgap, and renormalized bandgap are given by ∆, E 0 g , and E g . Note that this approach does not include cavity effects, which can be introduced in the Photon Green's functions solution following Ref. [19]. Quasi-periodic structures can also be addressed by a Green's functions formalism as shown in Ref. [28].

Numerical Application
The goal of this section is to illustrate the approach and the many quantities and parameters used making reference to published material, where the equations are used delivering very good agreement with experimental data. Photoluminescence is a very powerful tool to characterize semiconductor materials and map specific characteristics of new devices. Equation (58) is the reference, since it delivers the emission spectrum, which can be directly compared with experimental data. The carriers generated by the photo excitation process modify the spectrum and these modifications are described in Equation (58) approximately by the corrections induced by the (screened) potential, bandgap renormalization, and changes in linewidth governed by the dephasing or scattering Γ. The temperature T can be measured and used as an input parameter. The carrier density N can be estimated by measuring the input power and its spot size when focused on the sample, but in our recent investigations, where this theory has been very successfully compared with experiments [7][8][9], it has been treated as a free parameter, which has been globally adjusted. The other parameters that characterize the material, i.e., the fundamental band gap E 0 g , the electron and hole effective masses m e , m h , the static dielectric constant 0 , and the spin-orbit shift ∆, can either be found in the literature or robust numerical methods such as simulated annealing can be used to determine these parameters by direct comparison with experiments.
As a reference for the material parameters recently used for dilute nitrides and bismides and the corresponding bandstructure calculations that lead to the material parameters, see: Ref. [7] for GaAs 1−x Bi x ; Ref. [9] for InAs 1−x N x and Ref. [8] for more complex quaternary materials, such as InAs 1−x−y N x Sb y . The "s-shape" in the luminescence profiles as a function of temperature for these materials have been well explained in Refs. [7][8][9]. However, in the case of completely new materials, expecting to have general characteristics as the ternary or quaternary above, the parameters leading to the bandstructure may be unknown. This theory can be used as a numerical characterization tool as follows.
The corresponding bandstructure can depend on a number of unknown parameters for new compounds, but the approach used in Refs. [7][8][9] can be extended in the following way to extract these parameters by a systematic comparison between theory and experiments. For fixed excitation power, the luminescence can be measured for a number of different temperature points. The dephasing corresponding to different excitation processes can be calculated or taken also as a parameters. Thus, at each Temperature T, there is an ensemble of parameters, such as Experiments provide a series of data points measured at T = (T 1 , . . . , T N ). The calculated luminescence spectrum will be a function of T and will depend on the ensemble of parameters, denoted u C (T i ). The least squares method leads to estimates of the parameter ensemble E by minimizing the residual between the theoretical function and the experiments. Therefore, the problem becomes Trust Region-Reflective (TRR) methods deliver an efficient solution for this numerical problem and Ref. [8] gives further details of their application. Figure 5 depicts a numerical example to further illustrate choices for the main input parameters. Short period superlattices with strong delocalization of the electron and hole wavefunctions can be described in many cases by anisotropic 3D media, characterized by in-plane and transverse (along the growth direction) effective masses and dielectric constants. The anisotropy parameter γ is given by the ratio between the in-plane µ and perpendicular µ ⊥ reduced effective masses, γ = µ /µ ⊥ , with These can be calculated from the corresponding free carrier Hamiltonian, and full details of the method, which has led to good agreement with experimental data, can be found in Refs. [15,16,21,29]. Figure 5 shows calculated luminescence using the modified parameters determined by anisotropic medium theory for a short period GaAs-Al 0.3 Ga 0.7 As superlattice with repeated barrier and well widths equal to 2 nm. The resulting effective masses are m e ≈ 0.08; m h ≈ 0.12, m e⊥ ≈ 0.08, m h⊥ ≈ 0.53. These lead the anisotropy parameter γ = 0.67. The resulting exciton binding energy and Bohr radius are given respectively by e 0 = 5.37 meV and a 0 = 11.04 nm. Except of course for the actual value of the bandgap, which is larger due to quantum confinement, the strong delocalization of electrons and holes in this short period superlattice make the evolution of the luminescence with temperature look qualitatively similar to a three dimensional (bulk) semiconductor, notably the evolution of the line-shape. This is quite similar to the calculations presented in Ref. [9], which are in very good agreement with the experiments discussed in Ref. [13].

Conclusions
Photon Green's function techniques have been used to study different types of luminescence over the years, but the solutions are typically numerically intensive and not accessible for experimentalists or non-specialists. To bridge this gap, a simple analytical solution of the relevant Green's functions was necessary and the approach described here meets those needs. Notably, the evolution of luminescence with temperature has been successfully compared for GaAsBi [7] and InAs(N,Sb) [8] and InAsN [9] dilute semiconductors. Furthermore, the approach has been used to predict nonlinearities in short period superlattices treaded as anisotropic three-dimensional media [14][15][16]. However, details of the mathematical steps needed to achieve the final formulas used in these publications have not been previously presented and they are given in this review, complemented by numerical results demonstrating the range of validity of the main approximations used in the development of the approach. To complete the picture, an application section shows how to use the method as a numerical characterization machine, and the main parameters needed in typical simulations are illustrated with results for luminescence of short period superlattices. Screening of the Coulomb interaction between electrons and holes is discussed by means of the Hulthén potential and the steps provided are suitable as a guideline to the study of other interacting potentials of interest. They can be followed for further development of a suite of algorithms for efficient and easily programmable numerical characterization tools for a host of new bulk materials or superlattices that can described be as effective 3D media using anisotropic medium approximations. Except of course for the actual value of the bandgap, which is larger due to quantum confinement, the strong delocalization of electrons and holes in this short period superlattice make the evolution of the luminescence with temperature look qualitatively similar to a three dimensional (bulk) semiconductor, notably the evolution of the line-shape. This is quite similar to the calculations presented in Ref. [9], which are in very good agreement with the experiments discussed in Ref. [13].

Conclusions
Photon Green's function techniques have been used to study different types of luminescence over the years, but the solutions are typically numerically intensive and not accessible for experimentalists or non-specialists. To bridge this gap, a simple analytical solution of the relevant Green's functions was necessary and the approach described here meets those needs. Notably, the evolution of luminescence with temperature has been successfully compared for GaAsBi [7] and InAs(N,Sb) [8] and InAsN [9] dilute semiconductors. Furthermore, the approach has been used to predict nonlinearities in short period superlattices treaded as anisotropic three-dimensional media [14][15][16]. However, details of the mathematical steps needed to achieve the final formulas used in these publications have not been previously presented and they are given in this review, complemented by numerical results demonstrating the range of validity of the main approximations used in the development of the approach. To complete the picture, an application section shows how to use the method as a numerical characterization machine, and the main parameters needed in typical simulations are illustrated with results for luminescence of short period superlattices. Screening of the Coulomb interaction between electrons and holes is discussed by means of the Hulthén potential and the steps provided are suitable as a guideline to the study of other interacting potentials of interest. They can be followed for further development of a suite of algorithms for efficient and easily programmable numerical characterization tools for a host of new bulk materials or superlattices that can described be as effective 3D media using anisotropic medium approximations.