(image) (image) [ Previous ] [ Next ]

2 Heteroepitaxy and Selective Epitaxial Growth

The shape and quality of material regions is essential for device performance in modern semiconductor devices. The constant miniaturization increases the need for techniques which enable the fabrication of material regions with precision down to atomic levels. In order to obtain thin crystalline films with precise thickness and extraordinary quality, single atomic layers can be grown on top of crystalline substrate with a technique called epitaxial growth. The field of epitaxy encompasses all processes where single crystalline layers are deposited with well-defined orientations. This chapter gives an overview of these processes and discusses in detail how they are modeled in simulation frameworks.

2.1 Fundamentals of Epitaxy

The term epitaxy refers to the ordered deposition of well-defined crystal layers on top of a bulk crystal (substrate). The incoming atoms follow the structure of the crystalline substrate, which enables the emergence of thin films of structural high perfection down to the atomic arrangement. Epitaxial growth of thin films can be achieved with chemical and physical deposition techniques. Important examples for chemical vapor deposition (CVD) are atomic layer epitaxy and metalorganic chemical vapor deposition (MOCVD). A widely employed physical vapor deposition technique is molecular beam epitaxy.

Epitaxy proceeds on multiple time and length scales. Atomistic processes affect surface morphology and at the same time also the large-scale topography. On the microscopic scale epitaxial growth is driven by the movement of individual atoms, which takes place on microscopic time (ns) and lengths scales (Å). On macroscopic time (s) and length scales (µm) epitaxial growth results in the formation of crystal facets which are planes of extraordinarily high surface quality. Due to the high controllability and quality of the fabricated microstructures, epitaxy has received a lot of attention from the material and engineering sciences. Consequently, a plethora of epitaxy processes has been investigated. The broad field can be roughly categorized in homoepitaxy and heteroepitaxy. The former involves the ordered and crystal structure-preserving growth of the same material on top of a crystalline substrate. The latter refers to the growth of one kind of crystal on top of another crystal which does not necessarily have the exact same crystal structure [20].

Heteroepitaxial processes are well-known for the formation of precisely controllable thin films on planar substrates. This configuration is widely employed in the fields of optoelectronics and photonics (e.g., LEDs), where planar epitaxial films of precise thicknesses are required [36]. Furthermore, epitaxial growth typically exhibits a strong orientation dependence, i.e., the growth rates depend on the exposed surface of the crystalline substrate. The orientation dependence can be utilized to grow non-planar films which allow for the fabrication of complex 3D semiconductor geometries. Here, either a non-planar substrate or specialized masking of a planar substrate is required. If masking is employed, the epitaxy is desired to proceed selectively, i.e., a film is only deposited on top of a specific material, e.g., \( \silicon \) but not on the mask material, e.g., silicon oxide ( \( \siox \)). Many essentially equivalent terms for this mode of epitaxy are used in the literature: selective-area (epitaxial growth), selective growth, selective epitaxy, and selective epitaxial growth [37]. In this work, the latter term, abbreviated as SEG, is used while implicitly referring to non-planar SEG. A related growth mode is lateral epitaxial overgrowth, which is a special form of SEG with the property of the epitaxial film growing laterally over the mask [38].

In this section, epitaxy is addressed from several angles. First, the microscopic processes and crystal facets are discussed. Second, different growth conditions are introduced and finally, the growth techniques are presented.

2.1.1 Microscopic Processes During Epitaxial Growth

On the microscopic scale, the main processes governing epitaxial growth are adsorption, surface diffusion, nucleation, and film coalescence (illustrated in Fig. 2.1). The driving force of epitaxial growth is adsorption, where atoms or molecules originating from a gaseous or liquid atmosphere, referred to as adsorbate, adhere to a surface (substrate). The adsorbed atoms or molecules are generally referred to as adatoms or adspecies. The underlying mechanism of adsorption can either be physical (Van der Waals-type forces) or chemical in nature. The latter mechanism, called chemisorption, is based on a rearrangement of electron density and depends in general on crystal orientations and spatial distribution of adatoms. As such it is a kinetic process and can be quantified via the rate of adsorption \( R_\mathrm {ads} \). In a continuum picture \( R_\mathrm {ads} \) is proportional to the incident particle flux \( F \). The constant of proportionality is the sticking coefficient \( S \)

(2.1) \begin{equation} \label {eq:sticking} R_\mathrm {ads} = S \cdot F,\quad S\in \left [0,1\right ].                            \end{equation}

In general, sticking can be interpreted as the probability that an incoming particle reacts with the surface and adsorbs [39]. There is also a chance of desorption (i.e., the adatom is not incorporated into the crystal but enter the adsorbate phase again). On the microscopic level, the incorporation rate strongly depends on surface diffusion. A particle arriving at the substrate is not necessarily immediately incorporated into the crystal at a specific lattice site, but can move along the surface and hop from lattice site to lattice site. Thus, a redistribution of adatoms takes place via a diffusive process.

Surface diffusion is connected with a third microscopic process: Nucleation. Lattice sites on the free surface of the substrate are not energetically equal, i.e., there are surface sites that are more likely to adsorb an incoming particle. As a consequence, small clusters of adatoms are formed. Furthermore, adatom clusters represent an energetically favorable region for the subsequently incoming particles. Thus, adatom clusters act as nucleation sites which are the seed for further agglomeration. Since there are multiple nucleation sites distributed on the surface, they all grow in size.First independently, but then they form larger agglomerates until they coalesce to form a uniform film. [20], [36], [40].

(image)

Figure 2.1: Schematic illustration of microscopic phenomena during an epitaxy process. The crystalline substrate constitutes a step-terrace structure for the incoming adatoms. These are adsorbed at terraces and can move on the surface via surface diffusion. Clusters of several adatoms are formed, because nucleation sites are energetically favored. Depending of the growth condition, multiple individual islands or coalesced films are formed.
2.1.2 Crystal Facets and Miller Indices

On a macroscopic level, adsorption, surface diffusion, nucleation, and film coalescence result in crystal facets. These are flat surfaces of a crystalline material or film and are characterized by homogeneous orientation and defect-free atomic arrangement. A well-known example are facets of gemstones, which are cut along specific crystal planes to achieve their shiny appearance. In order to describe crystal facets, the Miller index notation is used [41]. Miller indices are based on a crystal lattice coordinate system with basis vectors \( \vec {a}_1 \), \( \vec {a}_2 \), \( \vec {a}_3 \). As shown in Fig. 2.2a, a facet is characterized by its intercepts \( x_1 \), \( y_1 \), \( z_1 \) with the lattice axes. Since the lattice is a regular tiling of space, the intercepts are integers. In the Miller notation of crystal planes

(2.2) \begin{equation} \label {eq:miller_plane} (h\,k\,l) = \left (\frac {s_1}{x_1}\, \frac {s_1}{y_1}\, \frac {s_1}{z_1}\right ),\quad s=\operatorname {lcm}(x_1, y_1, z_1),
\end{equation}

where \( \operatorname {lcm}(x_1, y_1, z_1) \) refers to the least common multiple of the axes intercepts, \( h \), \( k \), and \( l \) are the smallest possible integers that uniquely identify the crystal plane. If intercepts are negative, the respective indices are denoted with overlines, e.g., \( (h\,\overline {k}\,l) \). The set of crystallographically equivalent planes is denoted with curly brackets \( \{h\,k\, l\} \). A similar notation is used for crystal directions. A direction \( \vec {d} \) is expressed in terms of the lattice coordinate system

(2.3) \begin{equation} \label {eq:crystal_direction} \vec {d} = x_2 \,\vec {a}_1 + y_2 \,\vec {a}_2 + z_2 \,\vec {a}_3, \end{equation}

where \( x_2 \), \( y_2 \), and \( z_2 \) are integers. Again it is desired to denote directions with the smallest possible integers, which is achieved with the Miller indices for directions

(2.4) \begin{equation} \label {eq:miller_direction} \left [u\, v\, w \right ] = \left [ \frac {x_2}{s_2}\, \frac {y_2}{s_2}\, \frac {z_2}{s_2} \right ], \quad s=\operatorname
{gcd}(x_2, y_2, z_2), \end{equation}

where \( \operatorname {gcd}(x_2, y_2, z_2) \) denotes the greatest common divisor. The set of crystallographically equivalent directions is notated as \( \langle u\, v\, w \rangle \). An important special case should be noted for the cubic lattice: In a cubic lattice all basis vectors are perpendicular \( \vec {a}_1 \perp \vec {a}_2 \perp \vec {a}_3 \) and have equal length \( |\vec {a}_1| = |\vec {a}_2| = |\vec {a}_3| \). As a consequence, a plane with Miller indices \( (h\,k\,l) \) is perpendicular to the direction \( \left [h\,k\,l\right ] \).

(image)

(a)

(image)

(b)
Figure 2.2: Illustration of the Miller index notation. (a) Identification of a crystal plane according to (2.2) in a cubic lattice (characterized by three mutually perpendicular basis vectors \( \vec {a}_1 \), \( \vec {a}_2 \), and \( \vec {a}_3 \)). (b) In trigonal or hexagonal crystal systems, the Miller-Bravais convention with four indices is employed. Three (linearly dependent) basis vectors \( \vec {a}_1 \), \( \vec {a}_2 \), \( \vec {a}_3 \) reside in the so-called basal plane, representing three-fold rotation symmetry in a convenient way. The c-axis  \( \vec {c} \) is perpendicular to the basal plane. Three examples for crystal planes in Miller-Bravais notation are shown.

In the trigonal and hexagonal crystal system an extension of the index notation is used, referred to as Miller-Bravais indices. Four indices are used to indicate planes \( (h\, k\, t\, l) \) and directions \( \left [u\, v\, t\, w\right ] \) [20]. Fig. 2.2b illustrates that the lattice basis is composed of four vectors \( \vec {a}_1 \), \( \vec {a}_2 \), \( \vec {a}_3 \), \( \vec {c} \), which represents more naturally the underlying lattice. The basis is linear dependent, \( \vec {a}_3 = - (\vec {a}_1 + \vec {a}_2) \), thus index \( t \) depends on \( h \) and \( k \): \( t=-(h+k) \). Three planes with the corresponding Miller-Bravais indices are indicated in Fig. 2.2b.

Crystalline films can be fabricated with various epitaxial processes, which are discussed in the following sections.

2.1.3 Equilibrium Epitaxial Growth

If epitaxial growth takes place in thermal equilibrium, the incoming particles have sufficient time to move to the energetically most favorable lattice site. Adatom surface diffusion proceeds on a significantly faster time scale than the deposition process. This growth regime is called equilibrium epitaxial growth and can be understood with concepts from thermodynamics. The solid and liquid/gaseous phases are characterized by associated chemical potentials \( \mu \) and the driving force for crystallization close to thermal equilibrium is a difference \( \Delta \mu \) between the two phases. Furthermore, crystallization requires that the solid phase must be energetically more stable than the liquid or gas phase. In order to create an interface between solid and liquid/gas phase, energy has to be expended that depends on the surface energy density \( \gamma \) of each phase and interface [20]. Depending on the relation of the surface energies of substrate, interface, and film different thermodynamic growth modes can be observed [42], [43]:

  • • Frank-Van der Merwe growth: atoms are more strongly attracted to the substrate than to themselves. Film growth proceeds layer-by-layer.

  • • Volmer-Weber growth: atoms are more strongly attracted to each other that to the substrate. As a consequence, islands are formed, the film surface is not planar.

  • • Stranski-Krastanow growth: Frank-Van der Merve growth for first few deposited monolayer, then proceeds as Volmer-Weber growth.

The surface energy density of crystals depends in general on the crystal orientation  \( \vec {n} \) of the surface. In equilibrium growth conditions the macroscopic shape of a solid is determined by the energy minimization principle. The equilibrium shape  \( \mathcal {V} \) minimizes the surface term of Gibbs free energy [20]

(2.5) \begin{equation} \label {eq:ecs} G_\mathrm {surf} = \oint _{\partial \mathcal {V}} \gamma (\vec {n}) \mathrm {d}A \to \text {min}, \end{equation}

where \( \partial \mathcal {V} \) denotes the surface. For a liquid with isotropic  \( \gamma (\vec {n}) \) the equilibrium shape is a sphere (minimum surface area) and \( \gamma (\vec {n}) \) of a crystalline solid typically has multiple minima and maxima. These correspond to certain crystal facets which can be observed in the equilibrium shape. In order to determine the geometry of equilibrium shapes, Wulff’s construction can be used. The construction is based on Wulff’s theorem, which states that the length of a vector drawn normal to a equilibrium crystal face is proportional to its surface energy [44].

Aside from thermodynamic considerations, microscopic phenomena during equilibrium epitaxy have been studied in the literature [45]. On the microscopic level, epitaxial growth is associated with the formation of steps and terraces (step flow growth), as formalized by the Burton-Cabrera-Frank (BCF) theory [46]. Steps and terraces are schematically illustrated in Fig. 2.1. With this framework, observations like two-dimensional nucleation or the formation of screw dislocations and off-axis (vicinal) surfaces can be understood. An important concept is the Ehrlich-Schwöbel energy barrier. Microsteps on crystal surfaces can capture adatoms which are diffusing on a terrace with a certain probability. The capture probability depends on the direction from which adsorbed atoms approach the step. Adatoms originating from an upper terrace have to overcome a higher potential barrier, due to the lower number of nearest neighbors in the transition state. Therefore, there is asymmetry in the system which is quantified by the Ehrlich-Schwöbel energy barrier [47]. Several morphological properties can be explained with these approaches and thus form the basis of modeling approaches (see Section 2.2.2 and Section 2.2.3). Examples for these complex properties are meandering steps due to Bales-Zangwill instability [48], [49] and complex surface reconstruction due to dangling bond reduction [50].

2.1.4 Kinetic Epitaxial Growth

If particle flux is high and thus the deposition rate is high, the constantly incoming particles hinder the diffusion process. Adatoms are immediately (on-site) incorporated into the substrate, even though the lattice site is not necessarily the energetically most favorable. This growth regime is referred to as kinetically driven growth. From a thermodynamic perspective, growth proceeds far away from thermal equilibrium. One important consequence is the increase of elastic energy in form of elevated film strain, which ultimately limits film thickness [51]. The macroscopic crystal shape is governed by kinetic parameters, i.e., the growth rates of crystal facets. Since surface diffusion plays only a minor role, the density of adatoms on each facet reach stationary values. The growth rates are mostly determined by the process-induced deposition rate and thus by reactor-related quantities, like growth temperature or partial pressure of precursors [52]. The evolving facets can change considerably with changing reactor conditions, which is in contrast to the equilibrium growth governed by surface energy minimization. However, if reactor conditions are held constant, the local growth rates mainly depend on crystal orientation [53].

(image)

(a)
   

(image)

(b)
Figure 2.3: Illustration of the impact of the initial surface convexity on epitaxial growth. (a) The growth of convex initial geometries is dictated by the slowly growing crystallographic planes, and (b) the growth of concave surfaces is dominated by fast planes.
Reprinted with permission from Toifl et al., IEEE Access 8, 115406 [25] (2020). © 2020 Authors, licensed under the Creative Commons Attribution 4.0 License (http://creativecommons.org/licenses/by/4.0/).

In kinetic growth conditions the observed crystal shapes, referred to as kinetic crystal shapes [54], differ from the equilibrium crystal shapes governed by (2.5). In the limiting case of purely kinetic growth, surface diffusion is negligible and a stationary growth rate can be assigned to every surface plane. The associated growth rate function \( V(\vec {n}) \) is a stationary description of crystal growth (i.e., \( V \) is not a function of time). If \( \gamma (\vec {n}) \) is replaced with \( V(\vec {n}) \), Wulff’s construction can be conceptually applied to predict the crystal shape in kinetic growth conditions. While purely kinetic growth is an extreme case (with equilibrium growth on the other end of spectrum), this model works well for later stages of crystal facet evolution, where the crystal facets are large compared to typical surface diffusion lengths. However, initial stages of growth under out-of equilibrium conditions cannot be well-described, because surface diffusion and film stress play an important role in this case [55]. Nevertheless, a stationary growth rate \( V(\vec {n}) \) provides the possibility to describe and understand involved topographical consequences. Extrema of \( V(\vec {n}) \) directly manifest by the appearance of the corresponding facet in the kinetic crystal shape. If epitaxy proceeds on a convex surface, the plane with a relative minimum in growth rate determines the topography. On a concave surface the converse is the case [20], [56], as illustrated in Fig. 2.3.

(image)

Figure 2.4: Physical quantities contributing to the epitaxial growth rate in the continuum picture. While the incoming flux directly affects the growth rate, a more complex interaction between adatom incorporation \( \tau     \), surface diffusion, chemical potential, and film strain, results in a net contribution to the growth rate. This interaction can be modeled with methods from thermodynamics and continuum mechanics, as discussed in Section 2.2.5.

It is important to note that heteroepitaxy generally proceeds in between the extreme cases of equilibrium and stationary kinetic growth conditions. The net growth rate in a general heteroepitaxy process depends on multiple quantities which are summarized in Fig. 2.4. Deposition, adatom incorporation, surface diffusion, film strain, and chemical potential have inter-dependencies and thus mutually affect each other. In general, epitaxial growth does not perfectly follow the ideal evolution predicted by Wulff’s construction or the kinetic crystal shape. Nevertheless, the ideal growth conditions serve as valuable approximation. While a comprehensive and thus more involved description of general growth conditions can be achieved with methods from thermodynamics and continuum mechanics (Section 2.2.5), many specific epitaxy configurations are known to be primarily controlled by deposition kinetics, e.g., SEG of Si and \( \text {SiGe} \) [53]. In these cases, the stationary kinetic growth approximation is well-suited to capture the essential features of epitaxial film growths.

2.1.5 Selective Epitaxy of SiGe

In order to illustrate the most important characteristics of heteroepitaxial growth techniques, this section focuses on the specific case of SEG of \( \text {SiGe} \) in a CVD reactor. Growth is induced by subjecting the substrate to an atmosphere enriched with the species to be deposited. The desired atmosphere is achieved with precursors gases. In the case of \( \text {SiGe} \) SEG, typical precursors are silane ( \( \mathrm {SiH_4} \)), disilane ( \( \mathrm {Si_2H_6} \)), dichlorosilane ( \( \mathrm {H_2SiCl_2} \)), or germane ( \( \mathrm {GeH_4} \)). In general, the reason for selectivity is a high interface energy between masking material (e.g., \( \siox \)) and the epitaxial layer, causing the layer to avoid contact with the mask. Nevertheless, nucleated polycrystals are also formed to some degree on the mask. In order to remove these undesired polycrystal during the epitaxy process, an etch medium is introduced. Typical examples are hydrogen chloride ( \( \mathrm {HCl} \)) with high partial pressure [23] or chlorine gas ( \( \mathrm {Cl_2} \)), which remove nucleated \( \mathrm {Si_{1-x}Ge_x} \) on \( \siox \) with high efficiency [53], [57]. Deposition and etching can also be carried out in a cyclic process [52]. Tuning the relative gas flow of precursors allows to engineer the composition of the grown \( \mathrm {Si_{1-x}Ge_x} \) alloy. In-situ doping (i.e., doping introduced during growth) is also possible by adding dopant-specific precursors. Usually process temperatures are kept relatively low to assure absence of \( \mathrm {Si} \) or \( \text {SiGe} \) islands [58]. At temperatures in the range 575 °C to 900 °C, high quality growth with rates in order of 10 nm/min to 100 nm/min have been reported [23], [52], [53], [57].

Under the presented growth conditions, epitaxy proceeds near the stationary kinetic growth regime. The characteristic property of SEG is the pronounced anisotropic growth rate which depends on temperature, exposed surface area [37], doping precursor, and loading effects. Loading refers to variable deposition rates over the wafer caused by variations in surface topography or temperature [59]. In experiments, characteristic crystal facets \( \{1\,0\,0\} \), \( \{1\,1\,1\} \), and \( \{3\,1\,1\} \) have been observed. The typical ratio of associated growth rates is \( R_{3\,1\,1}/R_{1\,0\,0} \approx 0.5 \) and \( R_{1\,1\,1}/R_{1\,0\,0} \approx 0.25 \) at 750 °C [57].

In practice, it is important to consider the mechanical strain which builds up in selectively grown heteroepitaxial films. Since the crystal structure of substrate and grown film do not match (different crystal structure and/or lattice constants), mechanical strain is incorporated. For instance, in the plane of the hetero-interface, the strain \( \varepsilon _\parallel \) is given by [60]

(2.6) \begin{equation} \label {eq:strain} \varepsilon _\parallel = \frac {a_f-a_0}{a_0}, \end{equation}

where \( a_f \) refers to the lattice constant of the heteroepitaxially grown film and \( a_0 \) is the lattice constant of the substrate. If the strain exceeds a certain threshold, either elastic relaxation mechanisms (formation of 3D islands - Stranski-Krastanow growth, intermixing of deposited and substrate material), or plastic relaxation mechanisms (increased defect densities and threading dislocations) occur [51]. The consequence is a limited layer thickness which can be safely deposited. Even if the growth conditions are optimized, crystal defects are always present. In the context of SEG, advanced techniques to reduce the defect density have been investigated. One example is the technique of aspect ratio trapping, where threading dislocations can be eliminated by defect trapping in oxide sidewalls [52], [61].

2.2 Modeling Epitaxy

After discussing the fundamentals of epitaxy in the previous section, several approaches of modeling epitaxial growth are reviewed in the following. Even though this work is mainly focused on SEG, a more general perspective on modeling techniques which is not limited to SEG is given in order to understand the merits and limits of different approaches used to capture the intricate process of heteroepitaxy. Practical models need to combine a method to describe the geometric evolution and a technique to calculate growth rates based on a physical/empirical model. Growth rate models vary considerably in scope and fields of applicability, but can be roughly categorized in morphology- and topography-centered. Morphology-centered approaches target the intricate surface of epitaxial films and include lattice-gas models, BCF theory, and island dynamics models (IDM). Topography-centered approaches focus on the larger-scale geometry of emerging crystalline films and are thus well-suited for continuum models. In the following sections, an overview of these approaches is provided.

2.2.1 Atomistic Lattice-Gas Models

The atomistic description of epitaxial growth deals with adsorption, surface diffusion, nucleation, and film coalescence on the atomistic level. Since movement of individual adatoms is considered, these approaches are of particular interest in surface physics. A plethora of atomistic models have been devised to capture a wide range of morphological phenomena. They are based on a atomistic treatment of the interaction of substrate and atmosphere, hence they are referred to as atomistic lattice–gas models. Individual adatoms are tracked with their movement on discrete adsorption sites [62].

A well-suited algorithm to capture the evolution of such a discrete system is the lattice kinetic Monte Carlo Method (KMC) [63], [64]. It is a stochastic method which is based on a list of transitions and associated rates. The transitions abstract adsorption, desorption and diffusion kinetics. All transitions are assumed to happen instantaneously, independently and to be Markovian (memoryless). A common approach is to make use of the results of transition state theory and model temperature dependent transition frequencies \( \nu \) with an Arrhenius law

(2.7) \begin{equation} \label {eq:kmc_frequency} \nu = \nu _{0} \cdot \exp \left ( -\frac {E_\mathrm {a}}{k_B T}\right ), \end{equation}

with prefactor \( \nu _0 \) and activation energy \( E_\mathrm {a} \). The next transition and time increment is randomly chosen, e.g., with the direct Gillespie algorithm [65]. The crystal lattice is implemented by restricting particle positions and considering nearest neighbors to calculate binding energies (which eventually affect the transition frequency).

Specific studies on SEG on flat substrates have been reported [66]–[68]. Although the characteristic facets can be well-reproduced, the computational efficiency is problematic if desorption rates are high [67]. Nevertheless, lattice KMC methods provide the capability to study larger systems than first-principles approaches (e.g., molecular dynamics), because the abstraction to finite states and transitions reduces the computational complexity. However, the Arrhenius parameters of all transition processes have to be calibrated. This is possible either by expensive auxiliary density functional theory simulations or by fitting the results to an extensive experimental data set.

2.2.2 Burton-Cabrera-Frank Theory

The BCF theory [46] is a semi-atomistic description, in which continuum and microscopic quantities are prevalent. The main objective is to understand surface morphology features in equilibrium growth and thus BCF theory is particularly applied in the field of materials science. Crystal growth is viewed as a flow of steps, because the presence of steps is critical to allow for epitaxy growth. Adatoms diffuse on terraces towards steps and are only incorporated into the crystal if they reach a step before they desorb and leave the surface. The rate of advancement of a step is reduced to a diffusion problem on the surface, which is captured by a diffusion equation combined with the continuity equation for adsorbed molecules. To illustrate, the respective model for a one-dimensional problem with the y-coordinate along terraces and steps (nucleation on terraces is neglected) reads [40]

(2.8) \begin{equation} \label {eq:bcf} -D_\mathrm {S} \frac {\mathrm {d}^2 n_\mathrm {S}(y)}{\mathrm {d}y^2} = F - \frac {n_\mathrm {S}}{\tau _\mathrm {S}}.   \end{equation}

\( \mathrm {S} \) refers to a certain species, \( D_\mathrm {S} \) denotes the surface diffusion coefficient, \( n_\mathrm {S} \) is the number of adsorbed species per unit area, F is the flux of reactants, and \( \tau _\mathrm {S} \) is the mean resident time of the adsorbed species. Models in the spirit of (2.8) have been employed to qualitatively predict the occurrence of growth spirals on growing and grown surfaces due to screw dislocations [69]. The BCF theory is the theoretical basis for further models, as presented in the next section.

2.2.3 Island Dynamics Models

IDMs are an extension of the BCF theory. Crystal growth is viewed as a progression of island growth. The evolution of islands and steps is modeled as a free boundary model, where boundary movement is governed by mass conservation [31], [70]. Since the random formation of nucleation sites is modeled, the model is partly stochastic. The main idea is to employ continuum surface tracking methods to enable the simulation of the evolution of complex 2D and 3D epitaxial films. The most common surface tracking method employed for IDM is the level-set method, which enables surface tracking based on the continuous level-set function \( \phi (\vec {x}) \) and is discussed in depth in Section 4. Islands and steps are discretized by the integer number of atomic layers  \( k \). Island growth is described by a smooth evolution of the level-set function in the substrate plane. In vertical direction (i.e., growth direction) the surface can only propagate in discrete steps, thereby modeling terraces and steps. Atom diffusion on a specific terrace (constant  \( k \)) is considered as a continuous process and is thus described with a continuous diffusion equation. The time evolution of adatom density  \( \rho \) evolves through

(2.9) \begin{equation} \label {eq:idm_adatom} \frac {\partial \rho }{\partial t} = F + \vec {\nabla } \cdot (D \nabla \rho ) - 2 \frac {\mathrm {d}N}{\mathrm {d}t}.   \end{equation}

\( F \) refers to external particle flux (i.e., deposition) and \( D \) is the surface diffusion coefficient. Adatoms nucleate to form new islands at rate \( \frac {\mathrm {d}N}{\mathrm {d}t} \). The nucleation is modeled by randomly placing circular islands of size 2, which is assumed to be the minimum number of atoms required for a stable island [71]. The nucleation rate is defined using a mean-field approximation \( \frac {\mathrm {d}N}{\mathrm {d}t} \propto \langle \rho ^2(x)\rangle \), where the average is taken over all lattice sites. Boundary conditions for (2.9) consider the presence of Ehrlich-Schwöbel energy barriers. The actual surface is propagated depending on the velocity of island boundaries, which is driven by the gradient of adatom concentration at steps (satisfying mass conversation).

The same concept can also be employed with the phase-field method (see Section 2.2.5). In this case, a continuous approximation of the discrete height function of the growing film is used, resulting in blurred steps [72], [73]. Additionally, there have been further proposals for IDMs: A semi-implicit front-tracking method using parametric finite elements [74] and a boundary integral method [75]. In terms of scope, IDM studies are focused on island morphology on a (large scale) planar substrate. A example application is the investigation of instabilities during epitaxial growth of quantum dots [31].

2.2.4 Multiscale Approach

A further modeling approach is to consider film growth simultaneously on the micro- and macroscopic scale. The idea of tackling a scientific or engineering problem on different scales (usually length-scale) is referred to as multiscale modeling. Continuum models operating on the macroscale are combined with discrete or piecewise continuous models on the subscale [76]. Sun et al. [77] presented a combination of KMC on the microscale with an IDM on the macroscale to simulate epitaxial growth in two dimensions. On the macroscale, an IDM is solved with the level-set method. The missing microscopic data is estimated from a microscopic model, i.e., normal velocity and local adatom density. To achieve that, a set of microscale computational domains have to be constructed, initialized, and simulated with KMC. In order to translate the results to adatom density in the continuum macroscale picture, the heterogeneous multiscale framework [78], which is a general framework to couple macroscopic and microscopic model, is employed. Even though the coupling mechanism between the modeling scales add a significant layer of complexity, a reduction of the computational costs of KMC simulations is possible.

2.2.5 Continuum Models

In continuum models of heteroepitaxy the atomistic nature of matter is neglected and microscopic processes are abstracted with distributed quantities and material properties (e.g., surface diffusivity and adatom incorporation time). Consequently, the intricacies of surface morphology are not considered. Instead, the higher-level geometry on the macroscopic scale, i.e., topography, is at the center of attention. In particular, the evolving crystal facets are ideally flat and have no surface roughness. Fig. 2.5 shows modeling approaches, which correspond to the classification introduced in Section 2.1.4, thus targeting different growth conditions and areas of applications. The extreme cases of purely thermodynamics driven equilibrium growths and stationary kinetic crystal shape growth allow for making use of Wulff’s theory via surface energy functions \( \gamma (\vec {n}) \) or growth rate functions \( V(\vec {n}) \). The general case in between these extremes needs to be approached with comprehensive continuum modeling frameworks that enable incorporation thermodynamics, surface diffusion, adatom incorporation, and film stress [51], [55]. The concrete formulation of the computational models is always connected with the surface tracking method employed to follow the growth in three dimensions. The phase-field method is well-suited to integrate the general continuum model framework and the level-set method is applicable for stationary kinetic growth.

(image)

Figure 2.5: Summary of continuum modeling approaches and their associated growth regimes. Equilibrium growth is well-described with the Wulff shape, while stationary kinetic growth can be modeled with the related kinetic crystal shape. For any other (less extreme) growth regime, heteroepitaxy models based on Cahn and Taylor are well-suited [51], [79].
Cahn Taylor Phase-Field Models

The foundation for phase-field based continuum models is the classical Cahn-Hilliard model [80], which has been later extended by Cahn and Taylor [79]. In this work, the associated models are referred to as Cahn Taylor phase-field models. The crystal surface is modeled as diffuse-interface region, which captures the region of interaction between a gas and a solid phase. A phase-field function \( \varphi \) models the molar fraction of the solid component of the associated two-phase system. The boundary between a solid and a gas phase is \( \{\vec {x} | \varphi (\vec {x}) = \varepsilon \} \), where \( \varepsilon \) denotes the interface thickness [81], [82]. In recent years, continuum modeling of heteroepitaxy has been studied and refined by Bergamaschini et al. [51], [83]–[85].

The main assumptions are close-to-equilibrium heteroepitaxial growth at high temperature and/or low deposition fluxes. The chemical potential \( \mu \) is well-defined and growth is driven by a combination of direct deposition (external particle flux \( F \)) and equilibrium growth (proceeding due to surface diffusion and free energy minimization):

(2.10) \begin{equation} \label {eq:phase1} \frac {\partial \varphi }{\partial t} = \underbrace {F |\nabla \varphi |}_\text {Deposition} + \underbrace {\nabla \cdot \left [M \nabla \mu
\right ]}_\text {Surface diffusion} \end{equation}

Surface diffusion is driven by a gradient in \( \mu \) associated with the adatoms. \( M \) is the surface mobility which is in general temperature-activated and needs to be appropriately defined to restrict motion within the surface region only (i.e., surface movement is only possible near the interface). The functional form of \( \mu \) depends on elastic contributions, (anisotropic) interface energy, and heteroepitaxial intermixing [51]. One concrete formulation has been presented by Masullo et al. [83] in the context of CVD growth of \( \text {3C-SiC} \) (see also Section 5.3):

(2.11) \begin{equation} \label {eq:phase2} g(\varphi ) \cdot \mu = \underbrace {\frac {\delta G}{\delta \varphi }}_\text {Equilibrium} + \underbrace {\varepsilon \tau \frac {\partial
\varphi }{\partial t}}_\text {Adatom incorporation} \end{equation}

In this case, \( \mu \) is defined by the sum of the variational derivative of free energy \( G \) with respect to \( \varphi \) and a term that captures crystal facet dependent adatom incorporation. (2.10) and (2.11) form a coupled system of differential equations. \( g(\varphi ) \) is a function to improve the numerical convergence.

The merit of this approach is the physically intuitive integration of additional phenomena. For instance, continuum mechanics can be incorporated to allow for the consideration of elastic and plastic relaxation, which affect \( \mu \), and thus the growth evolution. Cahn Taylor phase-field models are well-suited for near equilibrium conditions, which enable the accurate prediction of the growth of nanostructures, e.g., nanowires or nanomembranes [55]. However, the numerical treatment of Cahn Taylor phase-field models is challenging: The finite element method has been utilized to solve the associated system of partial differential equations, which requires the careful tuning of stabilization functions, optimized adaptive meshes, and scalable algorithms to solve the full model for 3D configurations [51].

Stationary Kinetic Crystal Growth

A continuum modeling approach is also well-suited to describe epitaxy in stationary kinetic conditions. This growth regime is the main focus of this thesis. Here, the growth rate is determined by surface reaction kinetics. The net effect of incoming adatom flux and incorporation can be represented by a crystal orientation- dependent growth rate function \( V(\vec {n}) \). In these conditions, \( V(\vec {n}) \) provides a complete picture of the growth and can by directly correlated with the experimentally observable net growth rate of facets. From a measurement point of view, techniques like atomic force microscopy (AFM) or cross-sectional transmission electron microscopy (TEM) are well-applicable to characterize the growth rate of crystal facets. Therefore, an empirical \( V(\vec {n}) \), which provides the growth rate of every possible crystal orientation \( \vec {n} \) can be constructed. A function with this property can be directly employed in various general purpose surface tracking methods, e.g., cell-based [86] and explicit techniques [87]. A particularly convenient approach is the level-set method, where the time evolution of the surface topography is determined by the level-set equation [16] (see detailed discussion in Chapter 4):

(2.12) \begin{equation} \label {eq:levelsetequation_V} \frac {\partial \phi }{\partial t} + V(\vec {n}) \|\nabla \phi \| = 0.   \end{equation}

The topography is implicitly described via the level-set function \( \phi = \phi (x,y,z,t) \), which is a function of the spatial coordinates and time. All physical information about the evolution of the surface is contained in \( V(\vec {n}) \). Formally, the level-set method is related to the phase-field approach, with the most significant difference being the description of the gas-solid interface. In the level-set method the interface is considered abrupt, while a diffuse-interface region is used by phase-field approaches. From a computational perspective, the level-set description is very efficient, because no coupled system of differential equations needs be solved. The implications of employing a stationary kinetic crystal growth model \( V(\vec {n}) \) in the level-set framework is discussed in depth in Chapter 4. Even though stationary kinetic growth conditions are an extreme case, a growth rate function \( V(\vec {n}) \) can act as reasonable approximation if epitaxial growth is not significantly impacted by surface diffusion dynamics. This is the case if crystal facets length scales are large compared to the surface diffusion length. An important example is SEG, where typically large facets are formed.

A simple (semi-)empirical \( V(\vec {n}) \)-based approach is valuable for engineering applications: Topography plays an important role in semiconductor device design, especially when 3D device structures are considered. The growth anisotropy in combination with wafer misorientation or lithography mask misalignment might lead to formation of complex undesired geometric features which can impact subsequent process steps, e.g., ion implantation.

2.3 Summary

Table 2.1: Summary of modeling approaches and publications targeting heteroepitaxy and SEG.
Physics/Chemistry Surface Tracking Application Targets Dimension Publications
Atomistic (Lattice-Gas) Discrete Interface Surface and Island Morphology 2D [62]
Atomistic (KMC) SEG 3D on Flat Substrate [66]–[68]
IDM Level-Set Island Growth, Steps 2D [71]
Island Growth 3D on Flat Substrate [31]
Phase-Field Microstructure: Islands, Steps 2D [73]
Semi-Implicit Microstructure: Islands, Steps 2D [74]
Specialized Parametrization Stability Analysis of Steps 2D [75]
Multiscale Model; Micro: KMC; Macro: IDM Micro: Discrete Interface; Macro: Level-Set Enhance Performance of KMC 2D [77]
Continuum (Cahn Taylor) Phase-Field Near-Equilibrium Growth, Composite Substrates 3D [51], [83], [88]
Continuum (Stationary Kinetic) Level-Set Topography Simulation, Semiconductor Device Design SEG 3D [25]

In this chapter, the fundamentals of heteroepitaxy and selective epitaxial growth were discussed. Equilibrium and kinetic epitaxial growth conditions were introduced and important characteristics of a heteroepitaxy process were illustrated. Furthermore, approaches to model these processes were reviewed. Tab. 2.1 summarizes these modeling approaches and provides a selection of associated studies from the literature. Atomistic techniques (lattice-gas model) target the morphology of epitaxially grown films. Semi-atomistic approaches (BCF, IDM, and multiscale models) enable the simulation of certain microscopic features, e.g., island formation, on a larger scale. In contrast, continuum models target the topography of epitaxially grown films, which is of interest for the design of semiconductor devices. In particular, the capability of the level-set method to simulate stationary kinetic growth was identified.

In the next chapter, the subtractive anisotropic wet etching process is discussed and modeling techniques are reviewed. Anisotropic wet etching results in topographies that are similar to faceted heteroepitaxial structures.