Inhomogeneous charging and screening effects in semiconductor quantum dot arrays

The electronic properties of quantum dot (QD) arrays are strongly influenced by the Coulomb interaction of electrons within the dot as well as with those in neighbouring dots. In this paper, we investigate this behaviour taking into account screening by a free electron gas in the vicinity of the QDs. We find pronounced effects for standard capacitance[ndash]voltage (CV) measurements of QD structures embedded in a pn-diode. In particular, we show that the three-dimensional nature of the problem is crucial for devices with low dot-density, whereas the self-consistency between electron depletion in the bulk layer and dot occupation is important for high dot-densities. The Coulomb interaction between the dots induces a broadening of the peaks in the CV characteristic which is comparable with the effect of disordered QD arrays, where we considered realistic size and position fluctuations obtained by a kinetic Monte Carlo simulation.

Here ρ bulk (r) refers to the bulk charges, whose density is given by the quasi-thermal occupation of the conduction and valence band, determined by the differences between the respective quasi-Fermi levels and the electrostatic potential 3D (r). 0 and r (z) denote the absolute and relative permittivities, respectively, where z is the growth direction. See [24] for details. The charge of the QDs is modelled by the charge density where e > 0 is the elementary charge, αµ (r) the electron wave function of the level µ in dot α and p αµ ∈ {0, 1} denotes the actual occupation number of the levels. It is assumed that the occupations of the electron states are locally in thermal equilibrium with the conduction band where T is the temperature, E Fn (z QD ) the quasi-Fermi level of the electrons at the position z = z QD of the QD layer and E QD αµ the QD level energy including the effects of the self-consistent potential 3D (r). The discreteness of charge demands that the potential is determined by the actual discrete occupation numbers p αµ ∈ {0, 1} [25], which represent an ensemble following the thermal distribution function.
The self-consistent solution of the non-linear 3D partial differential equation (1) constitutes a formidable numerical task. Here this problem is circumvented by linearizing Poisson's equation around the solution 1D (z) of the 1D Poisson equation where the QD charge has been approximated as a sheet charge σ QD averaged over the crosssectional area A of the sample σ QD = −e A αµ p αµ (5) and χ QD is the characteristic function given by χ QD (z) = 1/ h for z QD − h/2 < z < z QD + h/2 and χ QD (z) = 0 otherwise. h is the height of the QDs located in the plane z = z QD . The bulk density ρ 1D bulk (z) is evaluated self-consistently with 1D (z) together with the quasi-Fermi levels E Fn (z). In [26], 1D (z) was used to determine the capacitance, and a significant level of broadening was needed to obtain agreement with the experiment. Now we linearize ρ bulk (r) around ρ 1D bulk (z) and equation (1) becomes Due to the linearity of the Poisson equation can be written as a sum of the 1D potentials 1D (z) andˆ 1D (z), given by and the 3D part QQ (r), determined by Here BQ = 1D +ˆ 1D can be conceived as a background potential for a single QD electron, excluding all other QD charges, while QQ describes the interaction between the QD charges. Equation (9) is solved by the Green's function G(r, r ) using standard techniques. Treating the interaction between the dots as well as with the background potential BQ as a first-order perturbation to the energy levels in the QDs, we find Here E intr αµ < 0 are the single particle level energies measured from the conduction band edge E c0 . Throughout this paper the values E intr αµ from k · p calculations [6] for InAs QDs are used, which explicitly give the size dependence: describes the level charging energy given by QQ due to the interaction of electrons occupying the QDs. The elements of the inverse capacitance matrix are calculated from the Green's function via where we approximately set ψ αµ (r) 2 = 1/V α within the QD volume V α , i.e. assuming a homogeneous spread of the wavefunctions over the QD volume. This approximation is justified for the interdot interaction (α = β), as the wave functions of different dots essentially do not overlap, and the spatial separation between different dots smears out inhomogeneities on the length scale of a single dot. The intradot interaction (α = β) is determined in a complex way by the real shape of the QDs, which modifies the wave function [6,8], and by both exchange [27] and correlation effects [7,9]. However, we do not expect any qualitative changes in our results. While typical experiments are performed for areas with millions of QDs, we have to restrict our calculations to a smaller cross-sectional area A, containing typically a 100 QDs. To mimic the larger size, we use periodic boundary conditions, i.e. we assume that the area is surrounded by eight equivalent areas with identical QD charge distributions. This is possible as long as the typical screening length, i.e. the extension of the Green's function in lateral directions, is shorter than the base length of the quadratic area A to avoid artificial double-counting.
In conclusion, the self-consistent, non-linear problem of calculating the QD energy levels with respect to the 3D potential distribution has now been reduced to the following iteration scheme, model A: (i) finding the self-consistent, 1D solution 1D , E Fn and E Fp of the Poisson equation (4) together with the current equations (see [24]) for a given σ QD ; (ii) solving the 3D linearized Poisson equation (9) to find the Green's function; (iii) calculation of the inverse capacitance matrix; (iv) statistical realization of a new QD occupation configuration to find the QD sheet charge σ QD by changing the occupation of a single QD state by one electron; (v) start again from step 1 with the new σ QD .
Note that the occupation of the QD levels according to (3) is coterminous to the stationary solution of a respective Master equation [28].
As the screening in equation (9) depends on 1D , and thereby on the average charging of the QDs σ QD , the capacitance matrix is changed whenever the occupation of a single QD state changes. Therefore, the iteration scheme of model A is highly time consuming. For small QD surface charges, one may neglect the influence of the QD occupation on the screening for the interdot interaction φ QQ and operate with fixed capacitance matrix elements. Then the iteration scheme reduces to model B: which was introduced in [24]. Finally, the depletion layer capacitance is given by the bias dependence of the vertical electric field in the depletion layer where U is the applied bias. For real systems, both the size distribution and the spatial arrangement of the QDs influence the behaviour, as they affect the charging energies E intr αµ and interdot capacitance matrix elements, respectively.

Results for ideal QD arrays
To study the behaviour, we consider a typical QD layer embedded in the n-doped side of a pnjunction with doping density 3 × 10 16 cm −3 [17]. The distance of the QD layer to the pn interface is set equal to 400 nm and we assume high doping in the p-doped layer. The QDs are assumed cylindrical with a diameter of 15 nm and a height of 3 nm. Arrays of 10 × 10 QDs are simulated and the temperature T = 70 K is used. At first we neglect any size or position fluctuations, thus all QDs have the same size and are periodically ordered.  Figure 1 shows the result for a small QD density of 10 10 cm −2 . At zero bias (see figure 2), the depletion region ends before the QD layer is reached and the QDs are occupied by up to six electrons each. With decreasing bias the depletion region increases and thus the capacitance decreases. Whenever the QD levels E QD αµ cross the quasi-Fermi level the states become empty, and the charge transfer is seen as a bump in the capacitance. Below −4 V, all QD states are empty. We find that the results of both models A and B agree below −2.6 V, where less than four electrons occupy the QDs. For higher bias, the increased sheet electron density reduces the bulk electron density in the vicinity of the dot layer, see figure 2 for 0 V. This causes a reduced screening for higher occupation, which is taken into account only in model A. Therefore, higher charging energies and lower occupancy are observed in comparison with model B, which neglects these effects.    The impact of screening on the capacitance matrix elements is shown in figure 3. It can be seen clearly that the interdot capacitance matrix elements drop significantly with increasing bias due to increasing screening effects.
The charging energy E char αµ given by equation (11) from our 3D simulation (labelled by 3D) is shown in figure 4. The results of models A and B differ slightly. For equal occupation, e.g. at −2.7 V, the charging energy is larger for model A due to reduced screening. On the other hand, the total charging energy is larger for model B at biases where model B gives an enhanced occupation. For comparison, we have also shown the quantity eφ 1D (z QD ) (labelled by 1D), which corresponds to the charging energy in a 1D simulation similar to [29]. One clearly sees that the correct treatment of the 3D nature is crucial for a good description of the device at the dot density of 10 10 cm −2 . In contrast, the self-consistent treatment of the screening is of minor importance here, so that the simplified model B from [24] is applicable.
At higher QD densities, the differences between models A and B increase as the capacitance matrix elements between nearest neighbours become much larger. This is shown for a dot density of 10 11 cm −2 in figure 5. Figure 6 shows that the total charging energies are much higher than at 10 10 cm −2 (figure 4). Therefore the dots can only be occupied by at most three (model A) or five (model B) electrons. The main part of the charging energy results from interdot interactions and thus the improved treatment of screening in model A is important. Note that the 1D result for the charging energy, eφ 1D (z QD ), does not differ too much from the 3D version for model A. Thus 3D effects are less significant at higher dot densities when compared with the effect of a self-consistent treatment of screening. Figure 6 shows that the charging energy varies by approximately 200 meV, while the lowest level of the QDs is charged with two electrons. On first sight this resembles a large variation in single-particle energies due to different dot sizes or compositions. Nevertheless, a perfect QD lattice was considered here and thus the effect is purely due to the change in the electrostatic environment as the neighbouring dots become more and more charged. Thus the charging effects can be easily mistaken for an inhomogeneous dot distribution.

Disordered QD arrays from KMC simulation of growth
To study the impact of disorder, we consider realistic QD arrays, obtained from a KMC simulation.
Here the self-organized growth of QDs in the Stranski-Krastanov mode is simulated taking into account self-consistently the heteroepitaxial strain field, as described in detail elsewhere [30]. The model uses an event-based continuous time KMC scheme, where diffusion and deposition of adatoms on a square lattice are the relevant processes. Starting from a perfect wetting layer, atoms are deposited at random on the surface and are allowed to move on the surface with a probability p according to Arrhenius' law Here ν 0 = 10 −13 s −1 is the attempt frequency, k B is Boltzmann's constant and T G = 700 K is the growth temperature. E S = 1.3 eV is the binding energy to the substrate and E n = (n − cn )E b with E b = 0.3 eV is the binding energy to the n nearest neighbours. The number of neighbours at the position of the adatom after the jump is denoted by n and the coupling with adjacent sites is described by the coupling constant c = 0.2. For adatoms moving vertically, e.g. jumping over an island edge, a Schwöbel barrier E SW = 0.1 eV is taken into account. All values have been carefully chosen such that good agreement with experiments is obtained. E str is due to the isotropic elastic strain field induced by the islands. It is calculated selfconsistently every 100th MC step and derived from elasticity theory using a Green's tensor formalism. Due to the lattice mismatch of the material system, line forces F i appear on the island edges S, and act as the source of strain. The elastic displacements u i (r) can be expressed by means of the static Green's tensor G ij (r, r ) (15) and the strain energy density is given by [31] with elastic constants λ = 2.9 × 10 11 erg cm −3 and µ = 3.5 × 10 11 erg cm −3 corresponding to averaged Lamé coefficients of InAs [30]. Since the screening length is of the order of the configuration area, in the CV calculations we consider a system where four identical island distributions are combined. Figure 8 shows that both samples have very similar CV characteristics, except that the features appear at slightly lower biases for sample II, exhibiting smaller single-particle-level energies E intr αµ due to the increased average size. Similar to the homogeneous case for a QD density of 10 11 cm −2 , the QDs are only occupied by three electrons in average at zero bias. The range of bias where the ground level of the QD is filled with two electrons is of the same order as for the case of a regular array. Thus, the irregular QD distribution is not manifested in the CV characteristics.
Note that there is a pronounced minimum in the capacitance around −2.9 V, where the average occupation is 1.5 electrons per QD. This feature was also found for a different random distribution of QDs at the same density, while the fluctuations between −2.5 and −1.5 V depend sensitively on the individual random realizations in our relatively small sample. The occupation distribution of the QDs given in figure 9 does not show any significant difference compared with different bias. Thus, we currently do not have an explanation for the origin of this feature.

Conclusion
We have investigated the influence of the 3D inhomogeneous structure on the charging effects in QD arrays taking into account screening of the interaction by bulk electrons in surrounding layers. As a general trend, we find that the 3D nature is very important for the charging energies for low QD densities. On the other hand, the self-consistent determination of the screening is crucial for high QD densities where the interdot Coulomb interaction dominates the charging. The increasing occupation of the QDs leads to a broadening of the charging peaks in the CV characteristics, which masks effects due to size or position fluctuations in realistic samples.