Phosphorus and Gallium Diffusion in Ge Sublayer of In0.01Ga0.99As/In0.56Ga0.44P/Ge Heterostructures

Abstract

This chapter gives a short review on dopant diffusion in germanium and specifies the underlying mechanisms of diffusion that involve the point defects. Box-shaped diffusion profiles are discussed that may be described as the phosphorus diffusion controlled by doubly ionized vacancies. In this mechanism, the diffusion coefficient depends on the electron concentration. The particulars of P and Ga diffusion profiles in the Ga-doped substrate of In0.01Ga0.99As/In0.56Ga0.44P/Ge heterostructures for multilayer solar cells are discussed. To calculate the diffusion coefficient, two methods were used: the Boltzmann-Matano (version of Sauer-Freise) and the coordinate-dependent diffusion analysis. It is established that coordinate-dependent diffusion analysis, which involves drift components together with diffusion components for diffusion profile description, is more suitable for description of the experimental profiles in such structures near p-n junction. A strong influence of intrinsic electric field on the dopant diffusivity was detected.

1. Introduction

Impurity diffusion in semiconductors is one of the main processes for electronic device manufacturing, but on the other side, it could badly influence a semiconductor structure in multistage high-temperature electronic device manufacturing processes. Dopants, as phosphorus, at diffusion temperatures are ionized; therefore they actively interact with ionized lattice defects creating charged complexes. These complexes are formed and destroyed in the diffusion process that leads to the appearance of generation and recombination components in a continuity equation that describes a diffusion process [1, 2].

Germanium is an important element to development of semiconductor theories and practice, and also it is a subject of many diffusion process researches. In this chapter, we focus on a narrow question: phosphorus diffusion in germanium, one of the main dopant of this material. Descriptions of diffusion processes were developing simultaneously with research of the crystalline and defect structure of this material and with improving of dislocation-free crystal growth technology together with development of measurement techniques and mathematical description of diffusion processes. That is why results that are 40 or 50 years old could be significantly different from contemporary ones. All these questions are under study and development. Progress in the first principal calculations together with the development of experimental techniques such as atomic force and scanning tunneling microscopy that allows to distinguish individual atoms and their lattice position will lead to the refinement of mechanisms and characteristics of diffusion processes. Our goal is to present the available data and knowledge about diffusion of phosphorus in germanium, possibly noting the problems and limitations of the representations used.

2. Phosphorus diffusion: first steps

Phosphorus, as a p-element of the group V of the periodic table, is a shallow donor impurity in germanium. The first works on phosphorus diffusion are about 1952–1954 years [3, 4, 5], and their review is in [2, 6].

It was previously mentioned that III and V group elements have a smaller diffusion coefficient than other groups of elements, and changes are mostly due to the frequency factor D₀. This was explained by their smaller ionic radius [5]. However, for elements of V group in germanium, this tendency was not confirmed (unlike that in silicon). Phosphorus, for example, having smaller ionic radius than any other V group element, has a smaller diffusion coefficient. For all shallow dopants (except of B), the activation energy is estimated as about 2.5 eV, and it slightly increased with decreasing diffusion coefficient in the range of As—Sb—P [5].

For a long time, constant diffusion coefficients were used for a fixed temperature [2, 3, 4, 5, 6]. These results were fairly expected, as in the absence of a reliable dopant profile measurement method, the diffusion coefficient was determined by p-n-junction depth; therefore it is in

$$D_P=1.2\cdot \text{exp}\left(-\frac{2.5}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}$$E1

[5] and taking into account the semiempirical Langmuir-Dushman formula:

$$D_P=2\cdot \text{exp}\left(-\frac{2.48}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}$$E2

At the same time, Ref. [5] already mentioned that high phosphorus concentration can lead to errors in calculations because of a tendency of this element to segregate. The surface concentration was not determined in the [5]. Another problem revealed in [5] was deviation of experimental values of p-n-junction depth in Sb diffusion (as the most studied dopant) from calculated dependence of p-n-junction depth on time ($d\sim \sqrt{t}$) at large time values. Therefore for estimation of the diffusion coefficient, a low diffusion time was used. Decrease of a penetration depth against expected one was attributed to diffusant evaporation in the diffusion process. These problems connected with the integral nature of a method of D coefficient determination.

In [7], the phosphorus profiles were determined using layered etching and sheet resistance measurements. Profiles of P in Ge that were made by vapor phase diffusion process were obtained for two surface phosphorus concentrations: less than and more than intrinsic carrier density n_i and at four diffusion temperatures—600, 650, 700, and 750°C. This allows to characterize temperature dependence of D. At low surface concentrations, the profile is described by Fick law, and diffusion coefficient is

$$D_1=330\cdot \text{exp}\left(-\frac{3.1}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}$$E3

At high surface concentration profiles which were extended, later [8] a name “box shaped” appears. For diffusivity calculations, authors applied Boltzmann-Matano method [1]. A dependence of the diffusion coefficient on the local phosphorus concentration was discussed. For the concentration-independent part, there was an expression obtained:

$$D_h=0.01\cdot \text{exp}\left(-\frac{2.1}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}$$E4

Experimental data did not fit well into Arrhenius curves, especially for data at high phosphorus concentrations. With the temperature increase, the diffusion activation energy also increased.

Similar results were obtained in [8]. SIMS method was used for concentration profile measurements. Phosphorus diffusion was carried out at temperature range 600–910°С. Surface concentration of phosphorus was higher than 10¹⁹ cm⁻³; therefore all samples were showing “box-shaped” profiles. Boltzmann-Matano method also was used for evaluating the concentration dependence of P diffusivity. The observed concentration dependence was approximately in agreement with results of [7]. The strong concentration dependence in D was attributed to dependence of D on Fermi level or due to strain effects caused by the difference in ionic radius of P in Ge.

In [8], temperature dependence of phosphorus diffusivity was found as

$$D_1=\left(0.009\pm 0.025\right)\cdot \text{exp}\left(-\frac{2.1\pm 0.2}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}$$E5

However the data of the paper allowed to derive another D(T):

$$D_1=1.21\cdot \text{exp}\left(-\frac{2.53\pm 0.2}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}$$E6

In the later works, a diffusion coefficient was called “intrinsic” for material, in which a dopant concentration n < n_i at the growth temperature, and it was called “extrinsic” when n > n_i.

In Figure 1, there is the dependence of D on phosphorus concentration from [5, 7, 8]. Data of [8] were calculated by Eq. (6).

Figure 1.

Diffusivity dependence on phosphorus concentration.

Integral values in [5] are noticeably higher than intrinsic D in [7, 8]; however, it does not exceed D in [7, 8] for high phosphorus concentrations.

Surprisingly, the experimental papers [7, 8] did not take into account extrinsic diffusion and dopant diffusion models, suggested in 1968 [9] and developed later [10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21]. Since vacancy in germanium is mostly acceptor with charge state up to −3, then positively charged phosphorus ion makes Coulomb-coupled pair with a charged vacancy. Diffusion of such pairs goes faster, and it was expected that it is in direct proportion to charged complex concentration.

3. Continuum theoretical calculations of dopant diffusion in semiconductors

The most detailed theory that describes dependence of dopant diffusivities on vacancy concentration in different charge states can be found in [10]. Indirect diffusion mechanisms, which involve vacancies V^k, are described by the following reaction:

$${\left(\mathit{PV}\right)}^j+\left(j-m-k\right)e^{-}=P^m+V^k$$E7

The local equilibrium is characterized by

$$\frac{C_{P^m}\cdot C_{V^k}}{C_{{\left(\mathit{PV}\right)}^j}\cdot n^{j-n-k}}=\mathit{\text{const}}$$E8

where P^m-phosphorus in substitution position, V^k-vacancies in k-ionization state, (PV)^j—vacancy-phosphorus complex in j-ionization state.

Generally, reaction (7) is a fast process compared to time scale of diffusion, which typically amounts to several minutes up to several hours. For this condition local equilibrium of the reaction is reached.

For the conditions near equilibrium:

$$D_P=\frac{C_{{\left(\mathit{PV}\right)}^j}{D}_{{\left(\mathit{PV}\right)}^j}}{C_{P^m}+C_{{\left(\mathit{PV}\right)}^j}}$$E9

If $n\approx C_P>n_i$,

$$D_{P_s^m}^{\mathit{eff}}=\left(m+1\right)D_{{\left(\mathit{PV}\right)}^j}{\left(C_{P_s^m}\right)}^{m-j}$$E10

Thus, for $m=+1,j=-1,D_P^{\mathit{eff}}\sim D_{{\mathit{PV}}^j}\cdot n^2$, if $j=-2,D_P^{\mathit{eff}}\sim D_{{\mathit{PV}}^j}\cdot n^3$

In one dimension, the diffusion equation takes the form:

$$\frac{\partial C_x}{\partial t}+\frac{\partial J_x}{\partial x}=G_x,$$E11

where C_x and J_x, respectively, are the concentration and flux of point defect X (P^m_s, V^k, P-V_j) as a function of time t and position x. Possible reactions between X and other defects are taken into account by G_x. If flux is determined by the diffusion of X, that is,

$$J_x=-D_x\frac{\partial C_x}{\partial x}$$E12

The diffusion equation is given by

$$\frac{\partial C_x}{\partial t}-\frac{\partial }{\partial x}\left(D_x\frac{\partial C_x}{\partial x}\right)=G_x$$E13

In [10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21], the behavior of P and Sb was consistently explained by means of the double ionized vacancy mechanism:

$${\mathit{PV}}^{-\mathit{1}}=P^{+\mathit{1}}+V^{\mathit{2}-}$$E14

If $m=+1,j=-1,\text{then}{D}_P^{\mathrm{eff}}\sim n^{\mathit{2}}$

In [17], As, Sb, and P were used for diffusion experiments. A Ge-dopant alloy source with about 1 at. % dopant content was used. Diffusion anneals were performed at temperatures between 600 and 920°C for various times in vacuum. The multiple use of the dopant source leads to depletion of the source. So the maximum doping level could be changed from the values that exceed the intrinsic carrier concentration n_i to values close or beneath n_i at the diffusion temperature. Doping profiles with penetration depths in the range of 30–150 μm were measured by spreading resistance method. Secondary ion mass spectrometry was used to record diffusion profiles with depths of a few microns. It was confirmed that in the range of low dopant concentration, the intrinsic diffusion with the constant D_in has been occurred. The extra diffusion with “box-shaped” diffusion profiles was observed when dopant concentration exceeded n_i. In this case:

$$D_{{\left(\mathit{PV}\right)}^{-}}^{\mathit{eff}}=D_{{\left(\mathit{PV}\right)}^{-}}\left(n_i\right){\left(\frac{n}{n_i}\right)}^2$$E15

$$D_p\left(n_i\right)=9{.1}_{-3.4}^{+5.3}\text{exp}\left(-\frac{\left(2.85\pm 0.04\right)\mathit{eV}}{k_{B}T}\right){\mathit{cm}}^{\mathit{2}}{s}^{-1}$$E16

Eq. (16) was calculated from Fickian-like profiles at low P concentrations. Then (15) were used for continuity equation, and a good agreement between experiment and calculations was achieved.

A “box-shaped” P profile was also detected under ion implantation procedure [18, 19, 20, 21]. The “quadratic model” was used to describe diffusion process.

In [21], the phosphorus distribution in germanium after ion implantation and annealing at temperatures 523 and 700°C was measured by SIMS method. It was shown that neither quadratic nor constant diffusion coefficient models cannot be used for profiles at 700°C annealing and longtime annealing for both temperatures.

Later a cubic dependence of the P diffusivity on the electron concentration was proposed [22]. The equations and dependencies used were.

$$\frac{\partial C_P}{\partial t}=-\frac{\partial J_P}{\partial x}{J}_P=-D^{\mathit{eff}}\cdot \frac{\partial C_P}{\partial x}-D^{\mathit{eff}}\cdot \frac{C_P}{n}\cdot \frac{\partial n}{\partial x}$$E17

$$D^{\mathit{eff}}=D^{\mathit{2}-}{\left(\frac{n}{n_i}\right)}^2+D^{\mathit{3}-}{\left(\frac{n}{n_i}\right)}^3$$E18

$$D_i=D^{\mathit{2}-}+D^{\mathit{3}-}$$E19

$$\begin{array}{l}{D}_i=44.3\cdot \text{exp}\left(-\frac{3.01\pm 0.04}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}\\ D^{2-}=11.1\cdot \text{exp}\left(-\frac{2.93\pm 0.01}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}\\ D^{3-}=5.7\cdot \text{exp}\left(-\frac{2.92\pm 0.02}{\mathit{kT}}\right){\mathrm{cm}}^2\cdot {\mathrm{s}}^{-1}\end{array}$$E20

There was a satisfactory conformity between experimental data and calculations for results of these authors and also with experimental data from [17] with this cubic model.

In Figure 2, a temperature dependence of the intrinsic diffusivity for cubic and quadratic models, experimental results in intrinsic diffusion regime [5] are presented. Figure 3 demonstrates concentration dependence D for two models together with experimental dependence [5] if proposed n = C_P. As we can see, calculated by Boltzmann-Matano values of D differ from estimations of D_in from Fickian’s part of diffusion curve, as it was done both in [17] for quadratic and [22] for cubic diffusion mechanisms.

Figure 2.

Intrinsic diffusivity for different models.

Figure 3.

Dependencies of diffusivity for cubic and quadratic models. Dashed lines are from experimental results [7].

4. Diffusion of phosphorus in InGaAs/InGaP/P heterostructures

In [16] co-diffusion of Ga and P was investigated, and it was shown that co-doping strongly affects the diffusion of phosphorus. The interest to Ga and P co-diffusion appeared with the developments of multicascade solar cells.

In last two decades, germanium is considered as the most suitable material for the first cascade of multiple solar cells based on A³B⁵ compounds that is intended for transformation of the infrared solar spectrum [23]. Germanium cascade of the multiple solar cells is formed by phosphorus diffusion into heavily gallium-doped germanium substrates. It was found that p-n junction depth weakly depends on the diffusion time. In [24, 25], P and Ga profiles in the heterostructure In_0.01Ga_0.99As/In_0.56Ga_0.44P/Ge were investigated. p-n junction of this element was formed at 635°C by phosphorus diffusion from In_0.56Ga₀.₄₄P buffer layer having thickness of about 24 nm to heavily doped of Ga germanium substrate (C_Ga = 2*10¹⁸ cm⁻³). The diffusion time was 2.6 min. SIMS has been applied to obtain profiles of P and Ga in heterostructure.

Figure 4 shows P, Ga, and free carrier concentration distribution in the Ge part of heterostructure. To calculate free electron concentration electroneutrality, equation was solved in the form of

$$C_P^{+}\left(x\right)+p\left(x\right)-n\left(x\right)-C_{\mathit{Ga}}^{-}\left(x\right)=0$$E21

Figure 4.

Profiles of P, Ga, n and p in Ge.

As dopant concentrations near interface are high, Fermi-Dirac distribution was used [26]:

$$n=N_C\cdot F_{\mathit{1}/\mathit{2}}\left(\eta \right),p=N_V\cdot F_{\mathit{1}/\mathit{2}}\left(-\eta -{\epsilon }_i\right)$$E22

where Fermi integral of order ½:

$$F_{1/2}\left(\eta \right)=\frac{2}{\sqrt{\pi }}\cdot \underset{0}{\overset{\infty }{\int }}\frac{{\epsilon }^{1/2}\mathit{d\epsilon }}{{\epsilon }^{\epsilon -\eta }+1};\epsilon =\frac{E-E_C}{\mathit{kT}};\eta =\frac{F-E_C}{\mathit{kT}};{\epsilon }_i=\frac{E_C-E_V}{\mathit{kT}}$$E23

where F is the Fermi level and E_c and E_v are bottom of the conduction band and top of the valence band, respectively.

Numerical calculations of Fermi level were made by Newton method for defined concentrations of P and Ga.

It was found that Ga diffuses insensitive to Ge substrates together with P. The higher solubility of Ga than P was found on the InGaAs/Ge interface as it was also noted earlier [27] that leads to formation of two p-n junctions. The shallow p-n junction was formed at a depth of 30 nm and the second one at a depth of 130 nm. Diffusion part of Ga profile demonstrated Fickian-shaped curve with D_Ga = 1.4 × 10⁻¹⁵ cm²/s that exceeds data 6 × 10⁻¹⁷ - 2.3 × 10⁻¹⁶ cm²/s [4]. As it was expected, phosphorus profile has two parts: Fickian type near the surface in p-region (C_Ga > C_P) and box-shaped between p-n junctions where n > n_i. Unfortunately using diffusion coefficient with quadratic and cubic dependencies, the P profile could not be accurately described [25].

Two methods of diffusivity calculations were used [28]. The first one was Sauer-Freise (SF) method based on the Boltzmann-Matano calculation of diffusivity [1]. The second one was method of the analysis of coordinate-dependent diffusion (CDD) [29].

In the CDD method, two parameters are introduced that describe a probability of hopping process ϕ(x) and probability that the nearest vacant place for diffusion is empty γ(x). Then diffusivity D(x) and drift velocity V(x) are expressed through these parameters and average distance between neighboring places λ. We have taken λ = a = 0.566 nm as a germanium lattice parameter. Then

$$D\left(X\right)=\varphi \left(x\right)\gamma \left(x\right){\lambda }^2$$E24

$$V\left(x\right)=\left(\gamma \left(x\right)\frac{\partial \varphi \left(x\right)}{\partial x}-\varphi \left(x\right)\frac{\partial \gamma \left(x\right)}{\partial x}\right){\gamma }^2$$E25

Drift term includes continuity equation:

$$\frac{\partial C_X}{\partial t}-\frac{\partial }{\partial x}\left(D_X\cdot \frac{\partial C_X}{\partial x}-V\left(x\right)C_x\right)=0$$E26

Figure 5 shows calculated dependencies of P diffusivity on x for both methods. Positions of p-n junctions are presented. As we can see, diffusivity calculated using SF method is comparatively higher than using CDD method. That may be a consequence of existing a strong electric field in the sample in the p-n junction regions that leads to appearance of a strong drift component in the charged particle diffusion.

Figure 5.

Diffusivity dependence on depth for T = 635°C.

Both methods of diffusivity calculations show two parts of D on x dependence: when x = 0–100 nm, diffusivity increases, and at higher values of x, diffusivity decreases. Width of the p side depletion region of the shallow left p-n junction on the Figure 5 is of the order of 5–8 nm $(\left|C_{\mathit{Ga}}-C_P\right|<{10}^{19}{\mathrm{cm}}^{-3}$); both sides of right p-n junction in Figure 5 are of the order of 50–80 nm ($C_{\mathit{Ga}}-C_p<{10}^{17}{\mathrm{cm}}^{-3}$); therefore an intrinsic electric field exists in the area between p-n junction. Approximately in the middle of junctions, the electric field changes its direction. Near the surface the intrinsic electric field accelerates negatively charged particles; when x > 100 nm, it inhibits diffusion. Outside of depletion regions (x > 160 nm), drift component of diffusion is negligible and diffusivity calculated by both methods which are equal.

Figure 6 shows dependencies of P diffusivity on n for Sauer-Freise, coordinate-dependent diffusion calculations, and different diffusion data from the literature.

Figure 6.

Diffusivity dependence on electron concentration for T = 635°C. 1: 0 < x < 25 nm, 2: 25 < x < 33 nm, 3: 33 < x < 60 nm, 4: 60 < x < 100 nm, 5: x > 100 nm.

An expected increase of the diffusivity with the free electron concentration was observed in both methods. Diffusivity produced by CDD calculation has two regions. The first one belongs to intrinsic diffusion (n < n_i = 3.2 × 10¹⁸ cm⁻³ [24]). As we can see, the lowest values of this part are equal to intrinsic diffusivity [5, 17, 30]. The second one corresponds to the diffusivity in the n-side of the p-n junctions and is higher than predicted both cubic and quadratic diffusion mechanisms. But the highest values of D = 2 × 10⁻¹² cm²/s at the n = 7 × 10¹⁸ cm⁻³ well correspond to maximum values [7], calculated by Boltzmann-Matano method. These values are observed in the electric field region of p-n-p structure that is formed in the germanium near the interface. Diffusivity dramatically drops at the ends of this structure in the p-region that may be connected with the shape of intrinsic electric field that in the case of linear p-n junction depend on x quadratically and drops sharply in the end of the depletion region. We can assume that the electric field causes not only the appearance of a drift component in diffusion but also increases the diffusivity of P-V pairs $D_{{\left(\mathit{PV}\right)}^j}$.

There are two regions of weak dependence of D on n. The first one at n < 2 × 10¹⁸ cm⁻³ corresponds to intrinsic diffusivity and is quite expected. The second is observed at high n in the region where the electric field exists. To understand the weak dependence of D_P on n (3–4 on Figure 6), we shall consider the equations for P-V complexes forming. In Table 1, the equations and parameters k and j that lead to different dependencies of $D_{{\left(\mathit{PV}\right)}^j}$ on n (see (7)) are presented for two cases.

	CP = n			CP = const
	k	j		k	j
n^0	−2	0	${\left(\mathit{PV}\right)}^0+e^{-}=P^{+}+V^{-2}$	−1	0	${\left(\mathit{PV}\right)}^0=P^{+}+V^{-1}$
	−3	−1	${\left(\mathit{PV}\right)}^{-1}+e^{-}=P^{+}+V^{-3}$	−2	−1	${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^{-2}$
				−3	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^{-3}$
n^1	−1	0	${\left(\mathit{PV}\right)}^0=P^{+}+V$	0	0	${\left(\mathit{PV}\right)}^0=P^{+}+V^0+e^{-}$
	−2	−1	${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^{-2}$	−1	−1	${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^{-1}+e^{-}$
	−3	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^{-3}$	−2	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^{-2}+e^{-}$
n^2	0	0	${\left(\mathit{PV}\right)}^0=P^{+}+V^{{\mathit{0}}^K}+e^{-}$	0	−1	${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^{\mathit{0}}+2e^{-}$
	−1	−1	${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^{-\mathit{1}}+e^{-}$	−1	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^{-1}+2e^{-}$
	−2	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^{-2}+e^{-}$
n^3	0	−1	${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^0+\mathit{2}{e}^{-}$	0	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^0+3e^{-}$
	−1	−2	${\left(\mathit{PV}\right)}^{-2}=P^{+}+V^{-1}+\mathit{2}{e}^{-}$

Table 1.

Equations and parameters for different dependencies of $D_{{\left(\mathit{PV}\right)}^j}$ on n.

The first is the same as in [17] when n = C_P+; the second is for the case of C_P+ = const as it is in our samples between p-n junctions (see Figure 6). We propose that $D_{{\left(\mathit{PV}\right)}^j}\sim C_{{\left(\mathit{PV}\right)}^j}$ [17].

Assumption that n = C_P is valid in a material with one type of impurity. In a strongly compensated material, the concentration of free charge carriers is significantly lower than the concentration of the impurity. Between p-n junctions in the measured heterostructure, phosphorus concentration changes slowly and we may suggest C_P+ = const. Phosphorus atoms in substitution positions are fully ionized, so m = +1; vacancies may be single, double, and triple ionized, that is, k = 0, −1, −2, and − 3, VP pairs—single and double ionized (j = 0, −1, −2).

Which type of reaction will be realized depends on the position of the Fermi level of the material, which controls the ratio of the centers in different charge state. The greater the electron concentration, the greater the charge state of acceptors, that is, for the condition n = C_P+, the most probable dependence of diffusivity of the complex is proportional n or n². In our case P-V complex should be charged, that is, j = −1 and −2. For C_P+ = const a weak dependence of the diffusion coefficient on n is possible most likely for the reaction:

$${\left(\mathit{PV}\right)}^{-1}=P^{+}+V^{-2}$$E27

The ionization energies of different charge states must be known to estimate a charge of a defect. It is obvious that ionization energies of vacancies and vacancy-assisted complexities depend on the temperature, but there are no reliable data of that energies [15, 31, 32, 33, 34, 35, 36, 37]. In [36] it was shown that at equilibrium conditions, half occupancy of the doubly negatively charged state of the vacancy-group-V-impurity atom pairs occurs when the Fermi level is situated at the middle of the forbidden gap. In spite of large phosphorus concentrations, n in the case of our interest is comparatively small, Fermi level is near the middle of the forbidden gap, and we may suggest that the (27) is an achievement.

As the electron density increases, the charge state of the pair can change. In the depletion region of the first p-n junction together with sharp increase of the Fermi level, the amount and charge of the pairs can be changed drastically, leading to a sharp increase in D_P.

5. Conclusions

In spite of numerous P in Ge diffusivity investigations, there are some issues that remain unclarified. The first one is the discrepancies between intrinsic diffusivities, calculated from Fickean type of diffusion profile at low phosphorus concentrations and those calculated using Boltzmann-Matano method from diffusion profiles at high P concentration. If we agree with vacancy assistant diffusion model, it means that P introduction into Ge increases the total vacancy concentration.

The formation of a p-n junction for germanium cascade of multiple solar cells due to the diffusion of phosphorus from the buffer layer In_0.56Ga_0.44P of In_0.01Ga_0.99As/In_0.56Ga_0.44P/Ge heterostructure leads to co-diffusion of P and Ga. The process was held at 635°C for 2.6 min. Solubility of Ga in the InGaP/Ge interface is higher than of P that leads to formation of two p-n junctions. Co-doping by gallium strongly affects the diffusion of phosphorus in germanium. We propose that it occurs primarily due to the electric field of the forming p-n junctions. P-type region is formed in the thin Ge surface layer (30 nm of order) with the depletion region thickness of 8–10 nm. The electric field of this p-n junction is directed to the Ge surface and accelerates both negatively charged Ga in interstitial positions and vacancy-phosphorus pairs. That leads to comparatively high gallium diffusivity D_Ga = 1.4 × 10⁻¹⁵ cm²/s.

We can point out that in the case of Ga and P co-diffusion, calculations of diffusivity by Sauer-Freise and coordinate dependence diffusion methods give values an order of magnitude higher than the values, obtained for quadratic and cubic diffusion model for phosphorus diffusion. An electric field of a depletion region of p-n junctions leads to the appearance of drift components of phosphorus diffusion. At low electron concentrations in p-region near Ge surface in which there is no an electric field, phosphorus diffusivity increases with n from intrinsic diffusivity values, produced from Fickean-type profiles at low P concentration, to that one calculated by Boltzmann-Matano method for high P concentrations, while P concentration sharply decreases. We may suppose the vacancy concentration increasing as the concentration of Ga and P that occupied the vacancies decreased.

It can be assumed that the electric field causes not only the appearance of a drift component in diffusion but also increases the diffusivity of P-V pairs. The sharp diffusivity growth and drop are consistent with the electric field direction. In the first p-n junction, it is directed to the surface and accelerates negatively charged particles including Ga⁻ and (PV)⁻. In the second one, it is directed into the sample that leads to decrease of the D(PV)⁻.

For a correct description of the Ga and P co-diffusion, it is necessary to take into account both changes in the concentration of charged centers due to a change in the Fermi level position and the formation and decay of diffusing pairs. For this, in the continuity equation, it is necessary to take into account not only the drift component but also the generation-recombination terms corresponding to the formation and decomposition of the diffusing pairs.