1. Introduction

The rise of diamond wire sawing requires changes in the texturing process of multicrystalline silicon substrates [1]. There is a strong need for an applicable structure with good optical performance, meaning low reflectance and strong light trapping. Black silicon structures have the potential to replace the isotexture as standard texture on multicrystalline silicon substrates. Among the possible fabrication processes for black silicon are metal assisted etching [1–3] and reactive ion etching [3–5]. Common to the resulting structures are high aspect ratios and size dimensions in the range or below the wavelength of the incident light.

A very low reflectance can result from structures with size dimensions far below the wavelengths of interest. This is due to the gradually increasing refractive index and can be modeled via effective medium approaches [6] [7]. However, for silicon solar cells not only a high front side transmission is relevant but also light redistribution into larger angles via scattering. This increases the path length of light in the silicon bulk and enhances the absorption close to the band edge. Therefore, integrating structures with characteristic dimensions in the range of the wavelength is beneficial due to the induced scattering. A prominent example is the current record solar cell on multicrystalline material [5]. However, those structure dimensions make optical modeling complex. Effective medium approaches and ray optical approaches are not applicable [6] [7]. The accuracy of scalar diffraction theory degrades substantially if the feature sizes get smaller than 10 λ, the refractive index is high, the slope angles are large and angles of incidence are increased [8]. Black silicon structures feature all of these properties, which makes scalar diffraction theory not a suitable model for describing them.

Although rigorous wave optical methods are computationally demanding, the general applicability of the rigorous coupled wave analysis (RCWA) to black silicon structures was shown by Kroll et al. [9]. They combine a coherent modeling of the black silicon structures with incoherent modeling of the substrate and obtain a good agreement between simulations and optical measurements of structured silicon substrates. Bett et al. introduced a cone approximation for black silicon [10]. They also reach a convenient agreement between pure RCWA simulations and experimental results of hemispherical reflectance and transmittance. The wavelength range of the study was however limited to 1000 nm - 1100 nm due to long calculation times of the RCWA for larger unit cells (up to 2 µm).

In order to overcome the problem of too long calculation times, field stiching methods for structures extending over large areas have been proposed for 1D [11–13] as well as 2D arrays [14]. These approaches use rigorous wave optical methods but replace a large unit cell with several overlapping subcells that are modeled instead. The electromagnetic field results are subsequently recombined from the central parts of the subcells and the diffracted order amplitudes are calculated via Fourier transforms.

In this work an advanced field stitching method is proposed that uses a field averaging algorithm in the overlap regions of the subcells. This allows the use of even smaller subcells and requires fewer subcell calculations since the information from the overlap regions is also used, in contrast to the information from the frame regions in [11–14]. The described modeling sequence requires less computational resources and can therefore handle larger unit cells. This is crucial when simulating random structures like black silicon. In addition to hemispherical reflectance and transmittance for a fixed angle of incidence also their angular distribution and an incidence angle variation can be modeled. This information can be used to determine redistribution matrices, which describe the angular light redistribution for a fixed set of angles and are used as input to OPTOS [15,16] simulations of complete solar cells or module stacks.

2. Summary of the field stitching approach

The simulation of large unit cells with rigorous methods is computationally demanding. The approach of this work is to simulate smaller subcells and calculate diffraction efficiencies of the original unit cell from a reconstructed field. It is an extended field stitching method (FSM) and includes the following steps, which are described in more detail and applied to a black silicon surface in section 3.

3. Modeling of black silicon surfaces

3.1 Generation of the height profile

The most intuitive approach to generate the height profile is an AFM measurement. However, within this work an approximation of the black silicon structure by cones was used, as proposed by Bett et al. [10]. This approximation was shown in combination with RCWA simulations to agree sufficiently well with experimental data. It furthermore allows a more systematic variation of structure parameters in future studies. The exemplary unit cell, which is investigated to demonstrate the validity of the FSM approach is depicted in Fig. 1(c). In the 2.4 µm × 2.4 µm unit cell 130 cones were randomly positioned with the bottom of the cones in the z₀ = 0 plane. This corresponds to a cone density of ~21-23 per µm² [17]. The cone diameters are also chosen randomly between 200 nm and 400 nm and the cone heights between 600 nm and 1000 nm. These parameters were chosen in agreement with reference [17] and the experimental structure in Fig. 1(a), 1(b). The z-layer thickness used in the RCWA simulations was 12.5 nm. Note that the FSM approach in principle allows choosing larger unit cells as well (see section 4). However, to be able to compare the results with RCWA calculations of the complete unit cell, a size of 2.4 µm × 2.4 µm was chosen in this section.

 

Fig. 1 a,b) SEM images showing a tilted view of the investigated black silicon structure in different resolutions. c) Tilted view of the modeled 2.4 µm × 2.4 µm structure unit cell, created by placing 130 cones with a maximum height of 1 µm in random positions based on the randomized cone approximation by Bett et al. [17].

Download Full Size | PDF

3.2 Deconstruction into overlapping subcells

In order to reduce the computational resources of the RCWA calculation the original unit cell was deconstructed into 16 subcells. These have an area of 1.0 µm × 1.0 µm each and include overlap regions with a width of 0.4 µm with each neighboring subcell. Overlap regions are required to enable a smooth field reconstruction between neighboring subcells, which will be explained in more detail in section 3.4. Simulating smaller subcells would in principle be possible as well. There are however restrictions. For example the simulation of single cones as basic unit cell led to strong deviations from the pure RCWA result of the original unit cell. This is attributed to the interaction between neighboring or even overlapping cones which can only be accurately considered if multiple cones are part of the same subcell.

3.3 Rigorous transmission modeling

As an example that is especially relevant for silicon solar cell applications, normal light incidence from air onto the structured silicon absorber was calculated for the complete unit cell as well as for each of the subcells individually. All of these simulations were done using periodic boundary conditions. The applied RCWA code is ‘Reticolo 2D’ [19]. Figure 2(a) shows the transmitted field result for the complete unit cell calculation at the bottom of the structure (z₀ = 0). The y-component of the field amplitude was chosen due to the incident polarization of the electrical field in the same direction, which leads to the largest amplitude values for the y-component. For the FSM calculation individual subcells are modeled. Figure 2(b) shows two neighboring subcells with overlapping areas visualized by a grey overlay. Note that the overlapping squares in Fig. 2(a) indicate the subcell positions but the depicted field information in (a) and (b) comes from different calculations.

 

Fig. 2 a) Transmitted field amplitude of 2.4 µm × 2.4 µm unit cell, calculated at once with the RCWA. b) Depiction of two neighboring subcells which were modeled individually for the FSM. The grey overlay visualizes the overlapping areas. The position of those subcells in the complete unit cell is indicated by the overlapping squares in (a).

Download Full Size | PDF

3.4 Reconstruction of complete field

The field information of the different subcells can be used to reconstruct an approximative complete field for the original unit cell. This was achieved by placing the fields of the subcells without overlap in their original position in the unit cell as depicted in Fig. 3(a) for the amplitude and Fig. 3(d) for the phase. Since all subcells were modeled individually, it is not surprising to find sharp discontinuities at the edges to the neighboring subcells. As these discontinuities are unphysical and can lead to an overestimation of the scattering into large polar angles, the field was smoothened. This was achieved by averaging the field in the overlap regions, where one value is available for each of the overlapping subcells. Introducing a linear weighting in the averaging process is the simplest way in order to give the field values close to the center of each subcell more weight than values close to the subcell border. The result of this linearly weighted average leads to the amplitude and phase in Fig. 3(b) and Fig. 3(e) respectively. Note the smooth transition between the neighboring subcells and the close agreement to the field, which was calculated for the complete unit cell without deconstruction into subcells in Fig. 3(c) and Fig. 3(f).

 

Fig. 3 Reconstructed electrical field (a) and phase (d) without averaging in the overlapping region and b,e) with weighted averaging in the overlapping region. c,f) Results of the RCWA modeling of the complete unit cell at once.

Download Full Size | PDF

3.5 Determination of angular intensity distribution

After reconstructing the field of the complete unit cell, the far field is calculated by Fourier transform for each component individually. Then, the absolute value of all components is used to determine the AID. For the transmitted fields described above, the AID is shown in Fig. 4. Note that his complete procedure is done for TE and TM polarization of the incident light and an averaging of the results afterwards.

 

Fig. 4 Polar angle dependency of the AID for normally incident light with a wavelength of 800 nm transmitted through a black silicon surface. Compared is the RCWA result of the complete 2.4 µm × 2.4 µm unit cell with the FSM.

Download Full Size | PDF

There are deviations between the complete RCWA result and the FSM. However, in the context of pathlength enhancement in silicon solar cells these deviations are not large. The question whether they significantly influence the absorption is evaluated in section 4 by comparing simulated with experimental results.

4. Comparison to experiment

For introducing black silicon into the OPTOS formalism, redistribution matrices were calculated via the simulation procedure described above for the 5.0 µm × 5.0 µm unit cell which is depicted in Fig. 5(a). This involves determining the reflection and transmission AID for a redefined set of angle channels (for more details concerning the channel definition, see [16]). For light incidence from air, the calculation was restricted to normal incidence since this is the direction that is relevant for the comparison with experiment. Within the silicon substrate, many different light directions occur due to scattering and diffraction. Thus, for light incidence from the silicon side, 25 distinct polar angles were modeled, all with the same azimuth angle. As the cones are randomly positioned no significant dependence of the redistribution properties of different incident azimuth angles is expected. Therefore, the results of the 25 angle pairs were used to fill the complete redistribution matrices. Note that the maximal simulation time for this specific unit cell, one incidence angle and one wavelength could take several days or even weeks using pure RCWA on a multi-core simulation computer. With the described FSM it was reduced to maximum several hours or a few days, strongly depending on the investigated wavelength. Figure 5(b) shows the transmission redistribution properties of the investigated structure for an exemplary wavelength of 1100 nm and normal light incidence from air. The light is strongly scattered into different angles, many of them outside the escape cone of a planar silicon air interface. Note that for determining the solar cell properties for different angles of incidence or in the module stack, also the light incidence from the front side of the structure would have to be calculated for all polar angles.

 

Fig. 5 a) Top view of the 5 × 5 µm² unit cell used for calculating the redistribution matrices. b) Representation of the transmission redistribution matrix for normal light incidence from air and an exemplary wavelength of 1100 nm.

Download Full Size | PDF

Applying the redistribution matrices, which were calculated via FSM, in the OPTOS formalism enables determining optical properties of substrates with black silicon structure. In Fig. 6, simulation and measurement results of a 250 µm thick silicon wafer with black silicon front side and planar rear side are compared. In the short wavelength range up to 1000 nm the reflection is very low for both the simulation as well as in the experiment. Above 1000 nm wavelength, reflection is slightly underestimated by the simulation. This also leads to a somewhat larger absorptance compared to the measurement. The only large difference appears for a wavelength of 1200 nm. Here, the neglect of free carrier absorption in the simulation could play a role. Having in mind the approximations of the modeling procedure, the agreement with the experimental results is reasonably good. Note that a regular approximation of the black silicon structure would feature low structure periods and thus a small number or even no higher diffraction orders that can propagate within the solar cell. Such a completely different angular distribution would result in a much lower absorptance compared to the measured values and the FSM.

 

Fig. 6 Measured reflectance, absorptance and transmittance of a 250 µm thick experimental sample with black silicon front side structure in comparison to simulation results, which were obtained using the FSM in combination with the OPOTS formalism.

Download Full Size | PDF

5. Solar cell application

The OPTOS formalism allows investigating also solar cell configurations with low additional effort once the redistribution matrices are available. The calculated black silicon interface can therefore be combined with other rear side matrices than the silicon-air interface, e.g. taking into account a metal reflector. In order to estimate the potential of this specific black silicon structure for silicon solar cells, a silver rear side mirror and a bulk thickness of 180 µm were chosen. The absorption result is depicted in Fig. 7 and compared to other state of the art textures and a Lambertian scatterer.

 

Fig. 7 Comparison of solar cell absorption simulations with different front side textures and a planar silicon silver rear side mirror. The black silicon structure does not only show a very low reflectance but also extremely good light trapping properties, which are due to the strong scattering behavior.

Download Full Size | PDF

The planar reference with anti-reflective coating (80 nm SiNx [20]) shows the expected large reflectance and low infrared absorptance. Introducing an isotexture as current state of the art texture for multicrystalline silicon solar cells (modeled based on an adapted spherical caps approach [21] with 80 nm SiNx ARC) reduces the reflectance and strongly increases the light trapping performance. The black silicon structure investigated in this work shows both a lower reflectance over the complete spectrum as well as a stronger light trapping performance for infrared wavelength. The performance is even close to a Lambertian scatterer, which was modeled using the same hemispherical reflectance and transmittance values as for the black silicon interface but with a Lambertian light redistribution.

To investigate the potential of black silicon for solar cell applications in more detail, structure parameters such as size, the distance or the aspect ratio have to be analyzed systematically. Furthermore, the encapsulated case has to be considered since the relative performance of different structures can change strongly on module level compared to the bare cell.

6. Conclusion

Within this work a black silicon structure is modeled optically using a combination of RCWA, a field stitching method, Fourier transforms and the OPTOS formalism. After dividing a large unit cell into smaller subcells the FSM uses RCWA results of the electromagnetic field in those subcells and a weighted averaging mechanism to reconstruct the field of a larger unit cell. This allows simulating unit cells as large as 5.0 µm × 5.0 µm in reasonable computation times. The comparison with optical measurements showed a reasonable agreement, slightly underestimating the reflectance in the long wavelength range. OPTOS simulations of solar cells revealed a significantly higher absorptance compared to a state of the art isotexture in the short as well as in the long wavelength range. The latter is due to the strong scattering of the black silicon which leads to a solar cell performance close to a Lambertian scatterer.