Silicon solar cells: basics of simulation and modelling: Using the mathematical program Maple to simulate and model a silicon solar cell

(1)

basics of simulation and modelling

Using the mathematical program Maple to simulate and model a silicon solar cell

Kisel solceller:

Grunderna för simulering och modellering

Sebastian Ekhagen

Fakulteten för hälsa, natur- och teknikvetenskap Kandidatprogramet i Fysik

22,5 Högskolepoäng

Handledare: Marucs Berg Examinator: Markus Rinio Datum Juni 2017

(2)

Abstract

The main goal of this thesis was to simulate a solar cell with the symbolic manipulation tool Maple and discuss the strength and weaknesses of using Maple instead of the already known simulation program PC1D. This was done mainly by solving the three essential differential equations governing the current density and excess electron and hole densities in the solar cell. This could be done easily by using known simplifications especially the low injection assumption. However it was also a success without using this particular simplification but the solutions had to be achieved using a numerical method instead of direct methods. The results were confirmed by setting up the same solar cell with PC1D. The conclusion is that Maple gives the user increased freedom when setting up the solar cell, however PC1D is easier to use if this freedom is not needed. At the end of this thesis a brief introduction is also made on the possibility of using Maple with a tandem cell setup instead of single junction.

(3)

II. Introduction to solar cells 4

A. Solar energy and AM1.5 . . . 4

B. PN junction in solar cells . . . 5

C. Generation and Recombination . . . 6

D. Surface recombination velocity . . . 9

III. 3 essential eqations for solar cells 10 A. The Current Density Equation . . . 10

B. Poisson’s Equation . . . 11

C. The Continuity Equation . . . 12

D. Summary . . . 12

IV. The Shockley-Queisser limit 13 A. Energy received from the sun . . . 13

B. Ideal Recombination . . . 13

C. Maximum Efficiency . . . 14

V. Solving the solar cell 16 A. Solving the depletion region . . . 16

B. Solving for the quasi-neutral regions . . . 19

C. Calculating the Efficiency . . . 20

VI. Solving for the qasi-neutral region with recombination simplification 22 A. Testing results with pc1d for low injection . . . 24

VII. Solving for qasi neutral regions without any recombination simplifications 29 A. Finding the efficiency under SHR-recombination . . . 34

VIII.From single junction cell to multi-junction cell 36 A. Tandem cell analysis . . . 38

IX. Concluding discussion 40 X. Appendix 41 A. List of symbols . . . 41

2

(4)

I. Introduction

Solar cells is one of many possible ways to harvest energy in a potentially environmentally friendly way. But to make it worthwhile and for countries to be able to make the transition to green energy the solar cells need to be cost effective and able to compete with fossil fuel plants. It’s important to notice that nearly all energy sources we use today on earth is produced from the sun, except of course nuclear energy and the Earth’s own thermal energy. For example oil is created by dead plants that once had grown by the use of solar energy. The reason that solar cells can be so efficient is because the time it takes for it to be renewed is no time at all, however the oil we use today was made 50 to 500 million years ago [1].

Solar cells uses the photovoltaic effect to produce energy from the sunlight and is thus called a photovoltaic device. This means that the device produces electric current when the material is exposed to light.

This is very powerful because there are no mechanical moving parts needed to produce power as is the case for windmills and hydro stations. Also no fuel is needed when in use so there is no environmental damage as a result.

Today’s problem with solar cells is to improve efficiency, the maximum efficiency for single junction solar cell is around 32% which will be calculated theoretically later on, but to even be able to reach this value of efficiency is nearly impossible and too expensive for it to be worthwhile for commercial use.

Thus solar cells used for residents in large scales are often in the efficiency range of 16% [2]. So it’s important to notice that every percentage of improvement of the solar cell’s efficiency is of huge value.

In labs the highest efficiency recorded is of 26.3%[3] but the methods used are then very expensive. So to study solar cells without having to produce them physically, programs are used to simulate possible architecture and to study what power output is possible with different input parameters. This is where programs such as PC1D are used, which will be used as a reference in this thesis. What PC1D does is let- ting researchers choose possible parameters such as the solar cells thickness, material, surface area and other useful inputs which will be covered in this thesis. This way the researcher can compare different setups of the solar cell to see which will be most cost-effective before having to produce them and thus saving money.

This thesis will discover the strength and weaknesses of using the symbolic manipulation tool Maple to simulate solar cells instead of using the already known PC1D and its competitors.

(5)

II. Introduction to solar cells

In this section we will study different phenomena and theories to be able to construct the physics of a solar cell.

A. Solar energy and AM1.5

The sun is always radiating photons with an energy of Ephoton = ~ω as a result of nuclear fusion of hydrogen in the core. For the sun (and many other objects) it can be approximated using black body radiation following from Planck’s law of black body radiation

I(ν,T ) = 2hν³ c²

1

e^hν/kT −1. (1)

The sun resembles very closely a black body with a temperature of 5800 degrees kelvin. But using this for solar cell physics would be too rough of an approximation due to the fact that photons from the sun have to penetrate our atmosphere. This results in a loss of some wavelength as they are absorbed by mostly water drops. So to make a better approximation specific radiation conditions are applied. The standard models are called AM0 AM1, and AM1.5, where the AM stands for "Air Mass". We will only be concerned with the AM1.5 which indicates that sunlight is radiating through the atmosphere at an angle of 41.8° above the horizon. The spectrum as a result of this assumption is shown in figure 1.

Figure 1:Spectrum of intensity Φ of photons from the sun hitting earths surface as a function of wavelength λ with the assumption AM1.5. Plot made in Maple using data from [12].

(6)

B. PN junction in solar cells

When two semiconductors, one of n doping (donor atoms are introduced) and the other of p doping (acceptor atoms are introduced) are brought together a PN junction is formed. This phenomenon is the whole basis of solar cells. In this section the phenomenon will be reviewed.

Again a PN-junction is formed when a semiconductor doped with donor atoms Nd and a semiconductor doped with acceptor atoms Na are brought into contact with each other, a boundary region between the n and p region is formed which we will call the space charge region. This region is formed when the holes of the p semiconductor and the electrons of the n semiconductor diffuse to opposite sides. However when this diffusion occurs the donor ions which are now left in the lattice will be positively charged as they have given away one electron, the same statement can be said for the acceptor atoms but these will be negatively charged. An electric field will now form in the space charge region, this electric field will work against the diffusion force such that the electrons on the n side will stay on the n side and same for the p side [5],[6],[7]. This is shown in the figure below.

Figure 2:Figure shows the PN-junction of a solar cell bombarded with photons.

The electrons will normally stay in the n-region due to the electric field which makes the electron experience a force F = −eE in the opposite of the diffusion, and the same happens for the holes which experience a force of F = eE in the p-region.It can also be seen in figure 2 that in this case the n semiconductor is hit by photons, this region is called the emitter the other is called the base. Usually the emitter is much narrower than the base. "It is more highly doped than the base by some orders of magnitude" [6].

It’s also important to notice that we will call the penetration depths of the pn junction is the place where the space charge region begins, xj. The width of the space charge region will be called W and the length of the solar cell will be called H.

(7)

C. Generation and Recombination

Generation of electron and hole pairs occurs when an energy source has given enough energy to an electron in the valence band such that it can "jump" from the valence band to the conduction band, where it can move freely across the crystal lattice. A solar cell uses mainly the photons from the sun as the energy source but can also receive energy from interactions with phonons as a result of thermal vibrations. The amount of energy needed to make the jump from the valence band to the conduction band is at the value of the energy gap of the material. For silicon the value of the energy gap is 1.124 eV which then means that a photon from the sun has to have a wavelength of ≈ 1.17µm maximum to be able to release an electron from the valence band [8]. We are only interested in electrons of the conduction band since they are the ones that can travel the circuit resulting in electricity.

However when we later study the real generation rate in a solar cell, we will take care of the re- flectance of the material and also the penetration depth of the photons in the given material. The generation rate is then calculated by the following formula:

G(x) = αN0(1 − Rf)e⁻^αx (2)

where α is the absorption coefficient of light in the semiconductor(α depends on the wavelength),N₀ is the photon flux at the face of the solar cell i.e the photon flux reaching the emitter shown in figure 2, which is calculated through N0 = Φ_hc^λ . Rf is the reflection coefficient which is found through:

R_f =

n₁−n₂ n₁+ n2

2 (3)

Here n1 is the refractive index of the first medium is always air, while n2 is the refractive index of the semiconductor. We can also see in eq.(2) the exponential decay of the absorption where x is the pene- trated length into the semiconductor.

Recombination however will work against the generation, which means that the electrons freely moving in the conduction band may lose their energy and recombine with the holes of the valence band while emitting a photon. This may occur due to several different phenomena which we will cover here.

1. Radiative recombination.

Radiative recombination is the most basic of recombination mechanisms. Radiative recombination is when an electron in the conduction band directly recombines with a hole in the valence band. In direct semiconductors the energy of the photon emitted from radiative recombination is at least the energy of the energy gap:

E_r = Eд= hc

λ (4)

However in indirect semiconductors the lowest possible energy (which is the gap between the lowest point of the conduction band and the highest point of the valence band) is increased due to the fact that the lowest energy point of the conduction band is not necessarily right above the highest energy point of the valence band. The transition between the conduction band and valence band in indirect semiconductor is more unlikely to happen because the electrons now has to interact with the crystal so the crystal momentum is conserved [5]. The recombination rate this way is directly proportional to the doping of the semiconductor which can be seen through the following:

(8)

R_r = Bnp (5) Where n and p are the concentration of free electrons and holes respectively and B is some proportionality constant. It’s stated in [5] that the proportionality constant is 10⁶ times larger than for the indirect band gaps. Thus proving that it’s the case that radiative recombination mainly occurs for direct band gaps.

2. Auger Recombination.

Auger Recombination is when an electron in the conduction band or valence band transfers its excess energy to an electron beside it. This makes the receiving electron excited whereby it loses its excess energy to collision with the lattice. This introduces phonons in the material which constitutes energy lost in the case of a solar cell. There are two ways the Auger recombinations occur which is an interaction between 3 particles, i.e hole-hole-electron in the valence band or electron-electron-hole in the conduction band [5]. For electron-electron-hole the recombination can be written as,

R_a = Bn²p (6)

and for hole-hole-electron

R_a = Bnp² (7)

Where B now is the Auger coefficient. The Auger coefficient B for silicon is 4 · 10⁻³¹ cm⁶s⁻¹[5].

3. Recombination due to trap levels.

When applying a dopant to the semiconductor one knows from basic solid state physics that energy traps can appear between the valence band and conduction band. For an n doped semiconductor the energy trap is close to the conduction band and for a p doped semiconductor this energy level is close to the valence band [9]. However when a dopant is introduced into the material it may give rise to defects, in the case of silicon this occurs when the electron configuration of the given atom doesn’t match that of a pentavalent or trivalent atom for n doping and p doping respectively [6]. These defects may then have an energy level which lie deeper in the band gap, this will get electrons and holes with energy of equal magnitude to be able to get stuck in these traps which then makes us lose possible conduction band electrons. This can easier be understood in the following picture which is also used in [6].

(9)

Figure 3:Different types of recombination due to trap levels. Holes are white circles and electrons are black. Et is the trap energy level. Adapted from [6].

We can see in figure 3 that there are four possible processes:

1. An electron is captured by a trap level from the conduction band.

2. An electron is emitted from the trap level into the conduction band.

3. A hole is emitted from the trap level into the valence band.

4. A hole is captured by the trap level from the valence band.

If every capture rate is calculated for each individual process, which is derived in [4], one finds with some calculation that the recombination rate of trap levels in one dimension is:

R_SHR = n(x)p(x) − n_i²

τ_p0(n(x) + n1)+ τn0(p(x) + p1) (8) Where τp0and τn0 is the mean carrier lifetime for holes and electrons respectively. Also n1and p1are

n₁ = niexp

− E_c −E_t k_BT

(9)

and

p₁ = niexp

− E_t −E_v k_BT

(10)

which are used to correspond to PC1D’s values found in [10]. Eq(8) is known as the Shockley-Read-Hall’s recombination rate.

(10)

D. Surface recombination velocity

As we discussed in previous section we saw that trap levels may occur in semiconductors due to defects as illustrated in figure 3. A solar cell have two surfaces: one in the front and one in the back, this results in increased defects due to the crystal lattice coming to an end. It follows that an increased trap level density at the surfaces results in an increased recombination rate, which in turn results in a shorter life time for the charge carriers. We then know that the excess minority carrier density lowers at the surface due to the increased recombination rate. We can thus describe this behaviour through the following expression evaluated at the surface:

−D_pd∆p

dx =S_p∆p (11)

From [5]. Where ∆p is the excess minority carrier density (in this case holes), Dp is the diffusion coefficient for holes. The above equation is important to be familiar with as this will be able to provide boundary conditions when we want to describe our solar cell during the simulation.

(11)

III. 3 essential eqations for solar cells

In this section the three essential differential equation that build up the solar cell will be derived.

A. The Current Density Equation

To understand the current density equation one needs to understand the meaning of drift and diffusion. In a solid electrons in the conduction band moves freely, but when no electric field is applied the movement of the electron is purely thermal and averages to zero in any direction [5].However if an electric field is present, electrons will be accelerated in the opposite direction of the applied field. This is shown in figure 4.

Figure 4: Figure is showing the drift of electrons with no electric field applied (black), and with an electric field applied (red). Adapted from [5].

The sudden changes in direction are due to electron-lattice interaction[6]. Also the movement of holes will be enhanced in the same direction as the applied electric field. With this in mind we can set up the equations for the current density due to drift of the electrons.

F = m^∗_ea =m_e^∗V_drift

t = −eE (12)

Where t is the time between collisions, e is the magnitude of the electric charge, m^∗_e is the effective mass of the electron and a is the acceleration. If we then take the meantime ¯t of all collisions we get:

V_drift = −1 2

eE¯t

m^∗_e (13)

We find the electron current density to be

J_ndrift= −neV_drift = ne1 2

¯t

m^∗_eE (14)

Equation (4) may be rewritten by introducing the electron mobility µ_e = 1

2 e¯t

m_e^∗ (15)

(12)

Which thus gives

J_ndrift = neµeE (16)

This is explained in more detail in [6]. By doing the same procedure for holes with a charge of +e we find the similar expression of the hole current density to be:

J_pdrift= peµhE (17)

Where

µ_h = 1 2

e¯t

m_h^∗ (18)

and m_h^∗ is the effective mass of holes.

However the drift is only one contributor of two to the current density in a semiconductor. The second contributor is the diffusion of electrons and holes. Diffusion for electrons and holes is the same as for example diffusion of air in a room, the molecules flow from a high concentration to a lower concentrations to reach equilibrium. The resulting current density from diffusion for electrons can thus be written as

J_ndiffusion = eDndn

dx (19)

and for holes

J_pdiffusion= −eDpdp

dx (20)

where Dn, D_p is the electron diffusion coefficient and hole diffusion coefficient respectively. The diffusion coefficient originates from kinetic theory and can be calculated by the Einstein formula which is

D = k_BT µ

e (21)

So with the drift and diffusion relations to current, the total current density for electrons can be written as:

J_n = J_ndrift+ J_ndiffusion = e(nµeE + Dndn

dx) (22)

and similarly for holes

J_p = J_pdrift+ J_pdiffusion = e(pµhE − D_pdp

dx) (23)

B. Poisson’s Equation

We know from basic electromagnetic field theory that free electric charges gives rise to an electric field through Poisson’s equation, this is also obviously the case for a semiconductor. We assume that the positively charged holes do not completely neutralize the negatively charged electrons. This in one dimension is thus

dE dx =

Q ϵϵ₀ = e

ϵϵ₀(p − n + N_d⁺−N_a⁻) (24)

Where p and n are the hole and electron concentration respectively, N_d⁺ and N_a⁻ are the donor atom concentration and acceptor atom concentration respectively. It’s important to understand that the donor

(13)

atom concentration counts as a positive charge since they give away an electron when introduced in the semiconductor, this makes the given atom a positive ion. The same argument follows for acceptor atoms as they produce a hole in the semiconductor, i.e takes up an atom making them negative ions.

C. The Continuity Equation

The continuity equation covers the contribution of generation and recombination of charge carriers. The equation in one dimension for electrons is written as:

d J_n

dx =−e(G(x) − R(x)) (25)

While for holes with a charge of +e is d J_p

dx =e(G(x) − R(x)) (26)

Where G(x) is the generation in one dimension and R(x) is the recombination in one direction. It’s intuitive to understand that the change of generations and recombination which is the value ^dn_dx for electrons and^dp_dx for holes. When these are multiplied with the charge of the charge carrier we of course find the change in current. For a further explanation consult [6].

D. Summary

The current density equation:

J_n = e(nµeE(x) + Dndn

dx) (27)

J_p = e(pµhE(x) − D_pdp

dx) (28)

Poisson’s equation:

dE dx =

e

ϵϵ₀(p(x) − n(x) + N_d⁺−N_a⁻) (29) The continuity equation:

d J_n

dx =e(R(x) − G(x)) (30)

d J_p

dx =e(G(x) − R(x)) (31)

(14)

IV. The Shockley-Queisser limit A. Energy received from the sun

The Shockley-Queisser limit is an important contribution to the field of solar cells which will be derived in this section. The Shockley-Queisser limit shows theoretically the maximum possible efficiency a single- junction solar cell can achieve under ideal circumstances.

Figure 1 will be used in this calculation. Assuming that the solar cell has its normal facing the sun, the number of photons at a certain energy gap hitting the solar cell per unit time, per square meter, per energy interval N=Eд can be calculated as:

N_=Eд = Φ E_photon ·

dλ dE =

Φ

E_photon · hc

E²_photon = Φhc

E³_photon (32)

From the calculation by [4]. Where Ephoton is the energy of the photon,Φ is found from figure 1. The total number of photons above a certain energy gap N>Eд reaching the solar cell is found through the following integral:

N_>Eд =

ż_4.42eV

Eдap

Φhc

E³_photondE_photon (33)

The maximum energy we can receive from photons from the sun is 4.42eV with the wavelength 280nm and will thus be used as our upper limit. It can be helpful plotting Eq(33) to get a picture of how the intensity of figure 1 relates to the number of photons hitting the surface above a given band gap energy.

Figure 5:The graph shows how many photons above a certain energy gap hit earths surface.

B. Ideal Recombination

Ideal recombination would be that the only possible way for electrons and holes to recombine is if they hit each other, this was named earlier as radiative recombination. When Shockley and Queisser first did the calculation they assumed that there were no Quasi Fermi Level splitting also known as QFL-splitting.

(15)

Without going into too much detail, it essentially means that when a bias is applied to the semiconductor the fermi-level for holes and electrons respectively split. No QFL-splitting is equivalent to the condition that the solar cell is under zero bias and in the dark. The recombination rate can then be calculated by the black body formula:

R₀ = 2π c²h³

ż∞ Eдap

E²_photon e^Ephoton^kBT −1

dE_photon (34)

Because a bias often is applied to the solar cell this recombination rate is not sufficient. When the QFL is splitting towards the conduction band by an energy E, the electron concentration is increased by exp(−E/kBt) (negative bias), while it is increased by exp(E/kBt) if the QFL splits towards the valence band (positive bias). In reality the splitting is more than the external applied voltage but for the ideal case a splitting of exp(eV /kBT ) [4], is assumed giving:

R_ra = 2π c²h³e^kBT^eV

ż∞ Eдap

E_photon² e^Ephoton^kBT −1

dE_photon (35)

This is though not entirely correct, this is explained in more detail by Würfel in [11], but it is sufficient for biases above 200meV.

C. Maximum Efficiency

The ideal current from electron-hole pairs is then easily calculated using the assumption form eq(35) that the only recombination occurring is from radiative recombination, and all possible photons are absorbed from eq(3) resulting in:

J_ideal(E_дap,V ) = e(N>Eд−R_ra) (36)

The maximum power can thus be derived by finding the maximum of the power Pideal, we derive Pideal

here:

P_ideal(E_дap)= V · Jideal(E_дap) (37)

then the maximum power Pmax satisfies:

d

dEP_ideal= 0 (38)

The maximum power was found numerically with Maple for each Eдap. The maximum efficiency is thus found through:

η_max(E_дap)= P_max(E_дap)

P_in (39)

P_in is the power per area reaching the solar cell, this value is found by integrating figure 1 this value is also called the solar constant. Which was found to be Pin = SolarConstant = 1000.36(W /m²). The famous Shockley-Quisser limit is found by plotting ηmax:

(16)

Figure 6: The famous Shockley-Queisser limit, showing the maximum efficiency that can be theoretically achieved at different band gaps.

By studying figure 6 we can see that for silicon solar cell that has a energy gap of 1.11eV can maximum produce an efficiency of η^Si_max ≈ 31, 8% while the maximum efficiency, of gallium arsenide at an energy gap of 1.43eV gets a maximum efficiency of η^GaAs_max ≈32.2%.

(17)

V. Solving the solar cell

In this section we will solve the three essential equations for a solar cell under different recombination conditions. This will be done by solving the equations in three parts of the solar cell which are, the depletion region and the quasi neutral regions i.e the emitter and the base regions.

A. Solving the depletion region

To start solving for the depletion region we will focus figure 2 on the said region and set up a new coordinate system for the given problem.

Figure 7:The figure shows the depletion region. Adapted from [3].

We can see in figure (7) that the width of the region of negative ions has the length xl and the width of the region of positive ions has the length xr. In figure (7) it seems that both regions have the same length, but this is only the case where the acceptor doping and the donor doping is of the same magnitude, but as we saw in previous sections this isn’t normally the case. Thus the length of both regions in figure (7) is dependent on the amount of doping atoms injected which will be derived as well.

If we take a look at Poisson’s equations once again:

dE dx =

e

ϵϵ₀(p(x) − n(x) + N_d⁺−N_a⁻) (40) This equation is not easily solved analytically stated by [6]. We thus need to introduce a set of assumptions to reduce the difficulty of the problem. These are given in [5] but will be rewritten here as well. The first major assumption is that the electric field is strong enough to drive any excess charge carriers out of the region thus n(x) = p(x) = 0. Also all dopants are ionized (Na = N_a⁻, N_d = N_d⁺), this is a perfectly fine assumption if we assume that the solar cell studied is at room temperature. We will also assume that negative and positive regions of figure (7) has abrupt edges, i.e no exposed ionized atoms are in the quasi-neutral regions xr ≤ x, x ≤ −x_l. The last assumption is that we assume that there is no electric field in the quasi-neutral regions due to that charge neutrality exist in these regions [5],[6],[8].

We can now directly see that eq.(29) difficulty dramatically decreases. Because n(x) = p(x) = 0 we get:

(18)

dE dx =

e

ϵϵ₀(N_d⁺−N_a⁻) (41)

Also the current density equations (17),(18) and the continuity equation (20),(21) becomes zero as a result of above assumptions. We choose to solve eq.(41) in the intervals −xl ≤x ≤ 0 and 0 < x < xr. Thus we get the equation for the interval −xl ≤ x ≤ 0 where N_d⁻ = 0:

dE

dx = − e

ϵϵ₀N_a⁻ (42)

and for the interval 0 < x < xr where N_a⁺ = 0:

dE dx =

e

ϵϵ₀N_d⁺ (43)

With the assumption that E(x) = 0 when x = −xl, x = xr gives us two boundary conditions to solve eq.(42) and eq.(43). We can now simply solve them as separable differential equations with the described boundary conditions.

E(x) = ż

− e

ϵϵ₀N_a⁻dx = − e

ϵϵ₀N_a⁺x + C1 −x_l ≤x ≤ 0 (44) And thus also

E(x) = ż e

ϵϵ₀N_d⁺dx = e

ϵϵ₀N_d⁺x + C2 0 ≤ x ≤ xr (45)

By evaluating the constants C1 and C2 with the assumption previously stated, E(x) = 0 when x =

−x_l, x = xr we get:

E(x) = −eN_a⁻(x + xl)

ϵϵ₀ −x_l ≤ x ≤ 0 (46)

E(x) = eN_d⁺(x − x_r)

ϵϵ₀ 0 ≤ x ≤ xr (47)

If we then plot the electric field with N_a⁻ = N_d⁺ = 5 · 10¹⁶ the relative permittivity of silicon ϵ = 11.9, xr

and xl are guessed to be 2µm in either direction, over the whole interval we get the following graph.

(19)

Figure 8:Electric field over the entire pn junction. (2µm=0.0002cm).

The electric field is continuous in the space charge region, this gives us, at x = 0:

N_a⁻x_l = N_d⁺x_r (48)

By studying above equation we see that the area of charges in the p region has to cover the same area of charges in the n region of the depletion region.

By now knowing the electric field through the whole depletion region, the potential can easily be found through

dφ

dx = −E(x) (49)

This will then be done in each region as eq.(42) and eq.(43) was done:

φ(x) = (−ş

CE(x)dx = şêNâ⁻_ϵϵ⁽^x−x₀ ^l⁾dx = êN_ϵϵâ₀⁻(^x₂² + xl ·x) + C3 −x_l ≤x ≤ 0

−ş

CE(x)dx = − ş^eN^d⁺_ϵϵ⁽^x−x₀ ^r⁾dx = −^eN_ϵϵ^d₀⁺(^x₂² + xr ·x) + C₄ 0 ≤ x ≤ xr (50) We are mainly interested in the potential difference of the depletion region which in turn lets us be able to set where the φ(x) = 0 is. So let’s say φ(−xl)= 0, we also know that just as for the electric field, the potential is continuous throughout the depletion region i.e φ(0) for −xl ≤ x ≤ 0 is equal to φ(0) for 0 ≤ x ≤ xr. This results in two boundary conditions that allows us to solve the constants in eq.(50), which by quick calculation in Maple gives us the potential in the depletion region to be:

φ(x) =











12 eN_a⁻(x+xl)²

ϵ ϵ₀ −x_l ≤x ≤ 0

12

e(N_d⁺x(x−2 xr)−N_a⁻x_l²)

ϵ ϵ₀ 0 ≤ x ≤ xr

(51)

The plot over the whole depletion region with the same constants and assumptions as for the electric field is thus:

(20)

Figure 9:Electrostatic potential over the depletion region in a pn junction

It can be seen in figure (9) that the potential difference is just the value of φ(xr). This leads us to the following:

φ_{di f f} = φ(xr)= 1 2

N_a⁻x_l²+ N_d⁺x_{r 2} e

ϵ ϵ₀ (52)

Where φdi f f = ^k^B_e^T ln

N_a⁻N_d⁺ n²_i

. We are now able to determine the length of xl and xr which in turn will give us the length of the depletion region. By using the assumption that N_a⁻x_l = N_d⁺x_r in above equation gives us:

x_l =

s2ϵϵ0

e

(φ_{di f f} −φ_app)N_d⁺

N_a⁻(N_a⁻+ N_d⁺) (53)

and

x_r = s 2ϵϵ0

e

(φ_{di f f} −φ_app)N_a⁻

N_d⁺(N_a⁻+ N_d⁺) (54)

In eq.(53) and (54), we have assumed that an external voltage is applied φapp. The total length of the depletion region is then found by

W = xl + xr =

s2ϵϵ0(φ_{di f f} −φ_app) e

s N_a⁻+ N_d⁺

N_a⁻ ·N_d⁺ (55)

B. Solving for the quasi-neutral regions

In previous section we studied the solution to the depletion region, however this is only one third of the way and solution for the emitter and the base will be discussed here. We’ve already made some assumptions from the previous section that have to apply here as well. For example we said that E(x) = 0 in the

(21)

regions x = xl and x = xr and thus also for the whole region of the emitter and the base. This makes the three essential equations simpler as written below:













J_n = eDndn

J_p = −eDdxpdp dE dx

dx = 0

d Jn

dx = e(R(x) − G(x))

d Jp

dx = e(G(x) − R(x))

(56)

We can easily see that this can be simplified further by differentiating the current density equation for both electron and holes.

d J_n

dx =eD_nd²n

dx² (57)

d J_p

dx =−eD_pd²p

dx² (58)

Then by just inserting the continuity equation for both cases which gives us two second order differential equations as follows:

d²n

dx² = R(x) − G(x)

D_n (59)

d²p

dx² = R(x) − G(x)

D_p (60)

By inserting the equation of generation eq.(43) and the recombination of eq.(8) we get for the p side solving for the minority carriers of course.

d²n dx² = 1

D_n

n(x)p(x) − n²_i

τ_p0(n(x) + n1)+ τn0(p(x) + p1) −αN₀(1 − Rf)e^−αx

(61) and thus also for the n side

d²p dx² = 1

D_p

n(x)p(x) − n²_i

τ_p0(n(x) + n1)+ τn0(p(x) + p1) −αN₀(1 − Rf)e⁻^αx

(62) When these two differential equations are solved one can simulate the solar cell.

C. Calculating the Efficiency

The total current in the solar cell is then calculated by finding the current density in each region, "we note that only those charge carriers generated either in the space charge region or at a distance of one diffusion length from the p-n junction contribute to current. Only this region of a solar cell is ‘active’"cited from [6].

J_T(λ) = JE(λ) + JDP(λ) + JB(λ) (63) Where JT is the total current density, JEis the current density genrated in the emitter, JDP is the current density generated in the depletion region and JB is the current density generated in the base. These can be calculated through the following

J_E(λ) = Jn(λ) x = xj + W

J_B(λ) = Jp(λ) x = xj (64)

(22)

The current density generated in the depletion region is calculated through J_DP = e

ż_x_j_+W

xj

G(x)dx = e ż_x_j_+W

xj

αN₀(1 − Rf)e⁻^αxdx = eN0(1 − Rf)e⁻^αx^j(1 − e⁻^αW) (65) For a closer explanation consult [6],[7].

The non ideal diode equation under illumination is given by

I = I0(e^eφapp^{N kBT} −1) − IL (66)

where I₀ is the dark saturation current, IL is the photocurrent from the solar cell and N is the ideality factor which takes a value between 1 and 2 depending on the recombination [7]. When N=1 gives the ideal ideality factor which can be derived by assuming low injection, and an infinite solar cell, this is done by Götzberger in [6]. By the use of eq.(66) we can derive the short circuit current (ISC) and the open circuit voltage (VOC). The open circuit voltage is found when I=0 giving:

0 = I₀(e^{N kBT}^eVOC −1) − IL→V_OC = Nk_BT e lnI_L

I₀ + 1

(67) and the short circuit current is found when V=0 giving:

I_SC = −IL (68)

The efficiency is then calculated through

η = F F J_SCV_OC

Φ_Light (69)

FF is called the fill factor and usually takes a value between 0.7 and 0.9 but can also be exactly calculated, this will not be done in this thesis but consult [6],[7] where FF is derived for a closer explanation.

(23)

VI. Solving for the qasi-neutral region with recombination simplification We can directly see that an exact solution to eq.(61),(62) are very difficult and direct application of the Maple solver does not produce any exact solution. Thus we need to try to simplify the equation further.

There is a possibility to simplify the recombination under certain assumptions and we will begin to use the assumption of low injection. Let us once again look at the Shockley Read Hall recombination.

R_SHR = n(x)p(x) − n_i²

τ_p0(n(x) + n1)+ τn0(p(x) + p1) (70) If we then assume that we are solving for the emitter first such that we have a high n doping. The holes are then our minority carriers p(x) we assume that n p for low injection. We will also assume that the total amount of holes is approximately the value of holes at equilibrium n ≈ n₀, and that n n₁ and p p1. Thus resulting in that eq.(70) now becomes for the n-side.

R_p = (n₀p(x) − n_i²)

τ_p0n₀ = (n₀p(x) − n₀p₀)

τ_p0n₀ = p(x) − p₀ τ_p0 = ∆p

τ_p0 (71)

By using the same arguments for the base where we have a p doping following from figure 2 where the electrons now are our minority carriers. We get the assumptions p n, p ≈ p0,n n1 and p p1

resulting in

R_n = (n(x)p₀−n²_i)

τ_n0p₀ = (n(x)p₀−n₀p₀)

τ_n0p₀ = n(x) − n₀ τ_n0 = ∆n

τ_n0 (72)

We can also assume that we are solving the equation for excess minority carriers ∆n, ∆p because the derivative of ∆n and ∆p are just n(x) and p(x) following from ∆n = n(x) − np0and ∆p = p(x) − pn0. If we now rewrite eq.(62) for the emitter and eq.(61) for the base we get:

d²∆p dx² = 1

D_p

∆p

τ_p0 −αN₀(1 − Rf)e⁻^αx

0 ≤ x ≤ xj (73)

d²∆n dx² = 1

D_n

∆n

τ_n0 −αN₀(1 − Rf)e⁻^αx

x_j + W ≤ x ≤ H (74)

The equations are now significantly easier and the solutions are found to be:

∆p = e⁻^Lp^x C + e^Lp^x D + α N₀ −1 + Rfe⁻^{α x}L_p² L_p²α²−1 Dp

0 ≤ x ≤ xj (75)

and similarly for the base

∆n = e^Ln^x B + e⁻^Ln^x A +α N₀ −1 + Rfe⁻^{α x}L_n² L_n²α²−1 Dn

x_j + W ≤ x ≤ H (76)

Where Lp = pDpτ_p and Ln = √

D_nτ_n, also known as the diffusion lengths which is the average length an electron in the conduction band diffuse until it recombines. To proceed we need to find two boundary conditions in respective region such that the constants A,B,C and D can be found.

We have already discussed boundary conditions at the surface due to surface recombination velocity following from eq.(11). We can also easily find one more boundary condition in respective region at the edges of the depletion region where we know from biased PN-junctions that p(xj) = p_n0e^eφapp^kBT and

(24)

n(x_j + W ) = np0e ^kBT . So to summaries the boundary condition for the setup following from figure 2.

For the emitter where holes are our minority carriers i.e (0 ≤ x ≤ xj):

−D_pd∆p

dx =S_p∆p(0)

∆p(x_j)= pn0e^eφapp^kBT −p_n0

(77)

Similarly for the base where electrons are our minority carriers in the region (xj + W ≤ x ≤ H):

∆n(x_j + W ) = n_p0e^eφapp^kBT −n_p0

−D_nd∆n

dx =S_n∆n(H) (78)

Below I will show the solutions to the differential equations for no applied bias as these are the only solutions that could in a reasonable way fit in the paper. However the general solutions will be plotted later for 0.2V.

By applying these boundary conditions we find that eq.(73) under no applied bias becomes in the emitter (0≤ x ≤ xj):

∆p = −αL²_p (−L_pαe^Lp^x −e⁻^αx)e⁻^Lp^xj + (e^Lp^xjL_pα + e⁻^αx^j)e⁻^Lp^xj −e^Lp^xje⁻^αx^je^Lp^xj (−1 + Rf)N₀

(L_pα + 1)((−Lpe^Lp^xjS_p + Dp)e⁻^Lp^xj + e^Lp^xj(S_pL_pe⁻^Lp^xj + Dp))(L_pα − 1) (79) and eq.(74) with Hb = xj + W becomes in the base (xj + W ≤ x ≤ H) :

∆n = (−1 + Rf) −αL_neH Lnα−Hb +x

Ln + αLneH Lnα+Hb +x

Ln + e⁻^{Lnαx+H −Hb}^Ln + e⁻^{Lnαx+H −Hb}^Ln −e⁻^{Hb Ln+H−Hb}^Ln −e⁻^{Hb Lnα+H−x}^Ln αN₀L²_n (−S_nL_ne^x−Hb^Ln + SnL_ne⁻^x+Hb^Ln + Dne⁻^{H +Hb}^Ln + Dne^{H −Hb}^Ln )(L²_nα²−1)

With these solutions we can now find the current density of both electrons and holes through the current(80) density equations thus resulting in the current density for holes in the emitter:

J_p = 4α(−1 + Rf)D_pL_pcosh_x

Lpj N₀ (L_pS_pe^Lp^x −D_p)e⁻^Lp^xj −e^Lp^xj(S_pL_pe⁻^Lp^x + Dp)₂

(L²_pα²−1)

·

− 1

2 (S_p(L_pα + 1)e^−αx + αDp)e^Lp^x −e^−αxD_pα L_pe⁻^Lp^xj + 1

2L_p(S_p(L_p −1)e⁻^αx−αD_p)e^Lp^xj + e⁻^αx^j(L_pS_pe^Lp^x − 1

2D_p)e^Lp^x + 1

2D_p(e^Lp^xje⁻^αxL_pα + e⁻^αx^je^Lp^x )

(81)

and the current density of electrons in the base:

(25)

J_n = D_n(−1 + Rf)αL²_nN₀ 2(L²nα²−1)L_nS_nsinh_H

b−x

Ln + Dncosh_H

b−H Ln

2·

·

e⁻H Lnα−Hb +x

Ln + e⁻H Lnα+Hb −x

Ln −αe⁻^{Lnαx+H −Hb}^Ln −αe⁻^{Lnαx+H −Hb}^Ln + 1

L_n(e⁻^{Hb Lnα+H−x}^Ln −e⁻^{Hb Lnα+H−x}^Ln ) + Sncosh(x − H_b

L_n )

−L_nα(e⁻H Lnα−Hb +x

Ln −e⁻H Lnα+Hb +x

Ln )+ e⁻⁻^{Lnαx+H −Hb}^Ln + e⁻^{Lnαx+H −Hb}^Ln −e⁻^{Hb Lnα+H−x}^Ln −e⁻^{Hb Lnα+H−x}^Ln

(82)

These are then our solutions for low injection, however to know if they are viable we need to do the same calculation in PC1D.

A. Testing results with pc1d for low injection

In this section we will discover if the results from previous section follows the same solutions as PC1D for silicon. We will first of all reduce the difficulty of the problem by making assumptions.

1. We will only consider one wavelength because if the problem can be solved for one wavelength then it can clearly be solved for any other wavelength under the same assumptions, at least at the level of solving differential equations. This assumption is made to remove any wavelength dependence from for example α and Rf in the generation eq.(2).

2. The voltage is applied through steady state.

3. PC1D takes into account several phenomena which are above the level of this thesis such as bandgap narrowing, electron and hole mobilities depends on doping density[10]. We will assume that the bandgap is constant and electron and hole mobilities are fixed values.

4. Constant temperature at 300K.

5.Uniform doping.

6.Monochromatic light of a wavelength λ=1140nm.

7. No contact resistances.

The constants used can be seen in the table below.

Constants

N_d 10¹⁹cm⁻³ T 300 K N_a 10¹⁶cm⁻³ Φ 0.1W/cm² n^si_i 1.08 · 10¹⁰cm⁻³ x_j 1µm µ^si_n 1110 cm²/Vs S_p 100 cm/s µ^si_p 410 cm²/Vs S_n 100 cm/s e 1.602 · 10⁻¹⁹C n₁ 1

ϵ^si 11.9 τ_n0 10µs

ϵ₀ 8.851 · 10⁻¹²F/m τ_p0 10µs k_B 1.381 · 10⁻²³J/K φ_app 0.2V

Table 1: Table showing general constants and material constants for silicon.

Where Sp and Sn is the surface velocity at the front and the back respectively, also τ_p0 corresponds to the emitter and τn0 corresponds to the base. Beyond the constants of table 1, the constants following

(26)

from our wavelength decision must be evaluated. These constants were taken from PC1D which had tables of them in the program[10]. We can see the following in the table below:

λ=1140nm

α=1 cm⁻¹ n₂=3.554 N0=5.73 · 10¹⁷ photons/cm²/sec Table 2: The table shows the wavelength dependent values.

Where n2 is the refractive index of silicon.

We can see in the figure below how the generation rate looks in the solar cell with the a monochromatic light of wavelength 1140nm. By using a wavelength with an energy below the band gap we can directly

Figure 10:Generation rate in a solar cell illuminated by a monochromatic light of wavelength 1140nm.

see in figure 10 that the generation rate in the solar cell is nearly constant due to the low absorption coefficient. We can now begin by plotting the solutions for the holes as excess minority charge carriers from eq.(79) resulting in the following:

(27)

Figure 11: Hole density for low injection assumptions in the emitter.

We can see from the figure above that the solution is only valid in the emitter region, that is (0 ≤ x ≤ x_j), which follows from the reasoning that this is indeed the only region where holes are the minority charge carriers.

We can then plot the solution for electrons through eq.(80) can be seen in the figure below. We can also

Figure 12: Electron density for low injection assumptions in the base.

see in this plot that its clearly understood that this solution is only valid for the region xj +W ≤ x ≤ H.

We can thus see in both figure 11 and 12 that the excess charge carriers goes quickly towards zero in the region close to the depletion region. Which of course follows from our boundary conditions that

∆p = pn0e^eφapp^kBT −p_n0 at x = xj and ∆n = np0e^eφapp^kBT −n_p0 at x = xj + W which are close to zero for a

(28)

low bias such as 0.2V. Let us also plot eq.(81) and (82) below: Notice that figure 13 is only plotted in the

Figure 13: Hole current density for low injection in the emitter.

region (0 ≤ x ≤ xj) as we saw in figure 11 this region is the only region where the current from minority charge carriers of holes makes sense.

Figure 14: Electron current density for low injection assumptions in the base.

Through the same discussion for the hole current density, we understand from figure 12 that we only need to look at the current density for the region (xj + W ≤ x ≤ H). Let us now use PC1D and tell it to solve for the same setup and plot the current density. The data received from PC1D was exported and

(29)

plotted in Maple.

Figure 15: The figure shows the comparison of the result received from PC1D (blue line) compared to the values I received (red line) from figure 14.

We can see in figure 15 that the results I received under low injection assumption corresponds very well the result received from PC1D. Only the electron current density was compared to one another as PC1D only chooses to plot the current density in the base. This means that the result of figure 13 cannot be compared. However this solution is only good for low injections.

(30)

VII. Solving for qasi neutral regions without any recombination simplifications In this section we will discover if we are able to solve the the 3 essential equations without simplifying the recombination as we did previously.

So let us once again write the Shockley Read Hall recombination rate R_SHR = n(x)p(x) − n_i²

τ_p0(n(x) + n₁)+ τ_n0(p(x) + p₁) (83) We can directly notice that if we want to solve for the whole recombination, the problem becomes significantly harder due to the n(x) and p(x) dependence. It’s thus no surprise that Maple couldn’t find an exact solution to the following coupled second order differential equations:











d²n

dx² = _D¹_n _n(x)p(x)−n2

τ_p0(n(x)+n1)+τn0(p(x)+pi 1) −αN₀(1 − Rf)e⁻^αx

d²p

dx² = _D¹_p _n(x)p(x)−n₂

τ_p0(n(x)+n1)+τn0(p(x)+pi 1) −αN₀(1 − Rf)e⁻^αx (84) To simplify the problem a bit we will assume no recombination at the surface and also no applied bias, this gives us for the n-side

D_pd∆p

dx = 0 for x = 0

∆p = 0 for x = xj

(85) and for the p-side

∆n = 0 for x = xj + W D_nd∆n

dx = 0 for x = H (86)

Again, we can’t find an exact solution to this problem. We thus need to solve it numerically, and the numerical method we will try to use is the shooting method.

The shooting method lets us use take our two boundary value problems and convert them into two initial value problems. This will be done by changing the boundary values at the depletion region edges to different values at the surfaces. So for the n region.

(∆p = 0 for x = xj) → (∆p = k for x = 0) (87) and for the p-side

(∆n = 0 for x = xj + W ) → (∆n = l for x = H) (88) We will then change the values of k and l such that the solution will "hit" our original boundary value at the depletion region edges. This is shown in the figure below for the n-region.

(31)

Figure 16:The figure shows an example on how the shooting method works

It can be easily seen in figure 16 that we will change the value of k in this case such that we find a solution to the boundary value problem of eq.(85). The same method will be used simultaneously for the p-side.

So let us begin setting up the constants such that we can see if we’re able to find a solution to the eq.(84), this time we will try to do so it closer relates to a real solar cell except the surface velocity condition of course. We will again use silicon as our material thus the material constants of table 1 will still be used.

The constants changed will be listed below:

Constants changed S_p 0 cm/s

S_n 0 cm/s φ_app 0V

Φ 0.092557 W/cm²

Table 3: Table shows constants change from table 1.

In addition to these constants, we will also change the wavelength and intensity of the wavelength such that it coincides with the spectrum of figure 1. We will use light of a wavelength of 750nm which holds an energy well above the band gap of silicon. The constants following from the chosen wavelength are listed below.

λ=750nm

α=1.3cm⁻¹·10³ n₂=3.723 N₀=3.50 · 10¹⁷ photons/cm²/sec

Table 4: The table shows constants following from a wavelength of 750nm.

As a wavelength of 750nm is used which carries a energy well above the band gap energy of silicion we can assume that the generation rate is no longer constant in the solar cell as was seen in figure 10.

The generation rate resulting from a wavelength of 750nm can be seen below.

(32)

Figure 17: Generation rate in a solar cell illuminated by a monochromatic light of wavelength 750m.

In relation to figure 10 we now see that the generation rate quickly goes to zero in figure 17 as a result of the absorption coefficient being large, which means that we only have generation in the beginning of our solar cell where all the photons are absorbed.

By then solving the coupled system of eq.(84) through the shooting method we find the following solution.

Figure 18:The figure shows the solution to excess electron density under 750nm monochromatic light in the base.

It is quite difficult to see if our solutions are valid in above figure. So let us zoom in at the x = xj+W