\chapter*{Acknowledgment}
\addcontentsline{toc}{chapter}{Acknowledgment}
-First of all, I would like to thank my official advisors Prof. Dr. Bernd Stritzker and Prof. Dr. Kai Nordlund for accepting me as a doctoral candidate at their chairs at the University of Augsburg and the University of Helsinki.
-I am gratefull to Prof. Dr. Bernd Stritzker who, although being an experimental scientist, gave me the opportunity to work on a rather theoretical field.
+First of all, I would like to thank my official advisers Prof. Dr. Bernd Stritzker and Prof. Dr. Kai Nordlund for accepting me as a doctoral candidate at their chairs at the University of Augsburg and the University of Helsinki.
+I am grateful to Prof. Dr. Bernd Stritzker who, although being an experimental scientist, gave me the opportunity to work on a rather theoretical field.
The working environment providing insight on problems dealt with in this thesis from the point of view of an experimentalist were of great importance contributing to the success of this work.
Furthermore, I would like to thank Prof. Dr. Kai Nordlund for his expertise in the field of atomistic simulations and the possibility to repeatedly visit his group in Helsinki.
Much progress of this study is owed to the regrettably few but fruitful discussions and explanations by electronic mail, within short meetings on conferences and, of course, during my long-term stays in the nice capital of Finland.
Next to living costs, travel expenses for scientific conferences as well as my stays in Helsinki were unbureaucratically taken over by the research foundation.
I owe my deepest gratitude to Prof. Dr. J\"org K. N. Lindner guiding me through my first time of scientific working already beginning with the diploma thesis.
-I am heartly thankful for his frequent encouragements and the confidence he has shown enabling a greatly independent research and the improvement of respective skills.
+I am heartily thankful for his frequent encouragements and the confidence he has shown enabling a greatly independent research and the improvement of respective skills.
Furthermore, his experiences in the materials system covered in this study and fine grasp with regard to scientific routines were of great help for this work.
Getting appointed to a professorship at the University of Paderborn, it was also him, who initiated a collaboration with the local theory group under the direction of Prof. Dr. Wolf Gero Schmidt.
At this point, I would like to express special thanks to Dr. Eva Rauls.
The present thesis would not have been possible without her assistance and mentoring with respect to the utilized methods required by the new approach of investigation as well as her expertise in the materials system.
I am greatly thankful for the possibility to repeatedly visit the theory group in Paderborn.
-In this context, Dr. Simone Sanna is acknowledged for respective technical support and Michael Weinl, doctoral student of Prof. J\"org K. N. Lindner back then, for accomodation.
+In this context, Dr. Simone Sanna is acknowledged for respective technical support and Michael Weinl, doctoral student of Prof. J\"org K. N. Lindner back then, for accommodation.
I am grateful to Priv.-Doz. Dr. habil. Volker Eyert for {\color{red}writing one of the certificates of this work.
-Furthermore,} his lectures on computational physics and the electronic structure of materials, which I attented during my academic studies, influenced me to pursue scientific research in the field of computational physics.
+Furthermore,} his lectures on computational physics and the electronic structure of materials, which I attended during my academic studies, influenced me to pursue scientific research in the field of computational physics.
-One more time, I would like to thank Prof. Dr. Bernd Stritzker for another two-month position as a member of his reaearch staff and various long-term employments as a reasearch assistant, which not only ensured a minimum of financial supply but also involved tutorships in the field of solid state physics that could be carried out in a more or less free and autonomous way.
+One more time, I would like to thank Prof. Dr. Bernd Stritzker for another two-month position as a member of his research staff and various long-term employments as a research assistant, which not only ensured a minimum of financial supply but also involved tutorships in the field of solid state physics that could be carried out in a more or less free and autonomous way.
-I am grateful to Ralf Utermann, responsible for the computing infrastructure in the physics department, for providing access and support with the excellently maintaind high performance units available in the Augsburg Linux Compute Cluster.
+I am grateful to Ralf Utermann, responsible for the computing infrastructure in the physics department, for providing access and support with the excellently maintained high performance units available in the Augsburg Linux Compute Cluster.
Furthermore, being employed as an assistant under his direction during the first times of my studies, he provided insight into modern computing technology and, in doing so, sparked my interest in computational physics.
I would like to show my gratitude to Rolf Anders.
-Similarly situated, performing numerical investigations at another experimental physics devision, he was always interested in discussing, both, specific and general problems related to physics as well as to scientific computing.
+Similarly situated, performing numerical investigations at another experimental physics division, he was always interested in discussing, both, specific and general problems related to physics as well as to scientific computing.
Residing in the office of the undergraduates, I would like to thank all of the many fellow students I met during this time.
Interesting and not necessarily scientific discussions appeared on the daily agenda.
It was a great time of lively exchange of ideas and excellent coffee that I will keep in good memory.
-I am likewise thankful to all the other members of the devision for the distinctive and pleasant working environment and the permanent willingness to discuss, which very often happened on fridays burning the midnight oil.
+I am likewise thankful to all the other members of the division for the distinctive and pleasant working environment and the permanent willingness to discuss, which very often happened on Fridays burning the midnight oil.
For discussions aside from natural science, I am grateful to the small but fine reading group formed on the fringes of the temporary student protests in 2009.
I would like to thank Meta Schnell, Berthold Arlt, Michael Lippok, Erika Rempel, Leo Sellinger, Michaela Strattner, Katja Teich and Matthias Link for the numerous sessions reviewing elements of the critical theory of society.
\noindent
Pierre Simon de Laplace phrased this vision in terms of a controlling, omniscient instance - the {\em Laplace demon} - which would be able to look into the future as well as into the past due to the deterministic nature of processes, governed by the solution of differential equations.
Although Laplace's vision is nowadays corrected by chaos theory and quantum mechanics, it expresses two main features of classical mechanics, the determinism of processes and time reversibility of the fundamental equations.
-This understanding may be regarded as the basic principle of molecular dynamics, considering an isolated system of particles, the behaviour of which is fully determined by the solution of the classical equations of motion.
+This understanding may be regarded as the basic principle of molecular dynamics, considering an isolated system of particles, the behavior of which is fully determined by the solution of the classical equations of motion.
\subsection{Introduction to molecular dynamics simulations}
-Molecular dynamics (MD) simulation is a technique to compute a system of particles, referred to as molecules, with their positions, volocities and forces among each other evolving in time.
+Molecular dynamics (MD) simulation is a technique to compute a system of particles, referred to as molecules, with their positions, velocities and forces among each other evolving in time.
The MD method was first introduced by Alder and Wainwright in 1957 \cite{alder57,alder59} to study the interactions of hard spheres.
-The basis of the approach are Newton's equations of motion to describe classicaly the many-body system.
+The basis of the approach are Newton's equations of motion to describe classically the many-body system.
MD is the numerical way of solving the $N$-body problem which cannot be solved analytically for $N>3$.
A potential is necessary to describe the interaction of the particles.
By MD, a complete description of the system in the sense of classical mechanics on the microscopic level is obtained.
\begin{enumerate}
\item A model for the interaction between system constituents is needed.
Interaction potentials and their accuracy for describing certain systems of elements will be outlined in section \ref{subsection:interact_pot}.
-\item An integrator is needed, which propagtes the particle positions and velocities from time $t$ to $t+\delta t$, realised by a finite difference scheme which moves trajectories discretely in time.
+\item An integrator is needed, which propagates the particle positions and velocities from time $t$ to $t+\delta t$, realized by a finite difference scheme which moves trajectories discretely in time.
This is explained in section \ref{subsection:integrate_algo}.
\item A statistical ensemble has to be chosen, which allows certain thermodynamic quantities to be controlled or to stay constant.
This is discussed in section \ref{subsection:statistical_ensembles}.
\end{enumerate}
-These ingredients will be outlined in the follwoing.
+These ingredients will be outlined in the following.
The discussion is restricted to methods employed within this study.
\subsection{Interaction potentials for silicon and carbon}
$U_1$ is a single particle potential describing external forces.
Examples of single particle potentials are the gravitational force or an electric field.
$U_2$ is a two body pair potential which only depends on the distance $r_{ij}$ between the two atoms $i$ and $j$.
-If not only pair potentials are considered, three body potentials $U_3$ or multi body potentials $U_n$ can be included.
+If not only pair potentials are considered, three body potentials $U_3$ or many body potentials $U_n$ can be included.
Usually these higher order terms are avoided since they are not easy to model and it is rather time consuming to evaluate potentials and forces originating from these many body terms.
Ordinary pair potentials have a close-packed structure like face-centered cubic (FCC) or hexagonal close-packed (HCP) as a ground state.
A pair potential is, thus, unable to describe properly elements with other structures than FCC or HCP.
-Silicon and carbon for instance, have a diamand and zincblende structure with four covalently bonded neighbors, which is far from a close-packed structure.
+Silicon and carbon for instance, have a diamond and zincblende structure with four covalently bonded neighbors, which is far from a close-packed structure.
A three body potential has to be included for these types of elements.
\subsubsection{The Tersoff bond order potential}
If the energy per bond decreases rapidly enough with increasing coordination the most stable structure will be the dimer.
In the other extreme, if the dependence is weak, the material system will end up in a close-packed structure in order to maximize the number of bonds and likewise minimize the cohesive energy.
This suggests the bond order to be a monotonously decreasing function with respect to coordination and the equilibrium coordination being determined by the balance of bond strength and number of bonds.
-Based on pseudopotential theory the bond order term $b_{ijk}$ limitting the attractive pair interaction is of the form $b_{ijk}\propto Z^{-\delta}$ where $Z$ is the coordination number and $\delta$ a constant \cite{abell85}, which is $\frac{1}{2}$ in the seond-moment approximation within the tight binding scheme \cite{horsfield96}.
+Based on pseudopotential theory, the bond order term $b_{ijk}$ limiting the attractive pair interaction is of the form $b_{ijk}\propto Z^{-\delta}$ where $Z$ is the coordination number and $\delta$ a constant \cite{abell85}, which is $\frac{1}{2}$ in the second-moment approximation within the tight binding scheme \cite{horsfield96}.
Tersoff incorporated the concept of bond order in a three-body potential formalism.
The interatomic potential is taken to have the form
\label{subsection:integrate_algo}
A numerical method to integrate Newton's equations of motion was presented by Verlet in 1967 \cite{verlet67}.
-The idea of the so-called Verlet and a variant, the velocity Verlet algorithm, which additionaly generates directly the velocities, is explained in the following.
+The idea of the so-called Verlet and a variant, the velocity Verlet algorithm, which additionally generates directly the velocities, is explained in the following.
Starting point is the Taylor series for the particle positions at time $t+\delta t$ and $t-\delta t$
\begin{equation}
\vec{r}_i(t+\delta t)=
\label{basics:verlet:taylor2}
\end{equation}
where $\vec{v}_i=\frac{d}{dt}\vec{r}_i$ are the velocities, $\vec{f}_i=m\frac{d}{dt^2}\vec{r}_i$ are the forces and $\vec{b}_i=\frac{d}{dt^3}\vec{r}_i$ are the third derivatives of the particle positions with respect to time.
-The Verlet algorithm is obtained by summarizing and substracting equations \eqref{basics:verlet:taylor1} and \eqref{basics:verlet:taylor2}
+The Verlet algorithm is obtained by summarizing and subtracting equations \eqref{basics:verlet:taylor1} and \eqref{basics:verlet:taylor2}
\begin{equation}
\vec{r}_i(t+\delta t)=
2\vec{r}_i(t)-\vec{r}_i(t-\delta t)+\frac{\delta t^2}{m_i}\vec{f}_i(t)+
\label{subsection:statistical_ensembles}
Using the above mentioned algorithms the most basic type of MD is realized by simply integrating the equations of motion of a fixed number of particles ($N$) in a closed volume $V$ realized by periodic boundary conditions (PBC).
-Providing a stable integration algorithm the total energy $E$, i.e. the kinetic and configurational energy of the paticles, is conserved.
+Providing a stable integration algorithm the total energy $E$, i.e. the kinetic and configurational energy of the particles, is conserved.
This is known as the $NVE$, or microcanonical ensemble, describing an isolated system composed of microstates, among which the number of particles, volume and energy are held constant.
However, the successful formation of SiC dictates precise control of temperature by external heating.
While the temperature of such a system is well defined, the energy is no longer conserved.
The microscopic states of a system, which is in thermal equilibrium with an external thermal heat bath, are represented by the $NVT$ ensemble.
-In the so-called canonical ensemble the temperature $T$ is related to the expactation value of the kinetic energy of the particles, i.e.
+In the so-called canonical ensemble the temperature $T$ is related to the expectation value of the kinetic energy of the particles, i.e.
\begin{equation}
\langle E_{\text{kin}}\rangle = \frac{3}{2}Nk_{\text{B}}T \text{, }
E_{\text{kin}}=\sum_i \frac{\vec{p}^2_i}{2m_i} \text{ .}
The advantage of the approach is that the coupling can be decreased to minimize the disturbance of the system and likewise be adjusted to suit the needs of a given application.
It provides a stable algorithm that allows smooth changes of the system to new values of temperature or pressure, which is ideal for the investigated problem.
-\section{Denstiy functional theory}
+\section{Density functional theory}
\label{section:dft}
-Dirac declared that chemistry has come to an end, its content being entirely contained in the powerul equation published by Schr\"odinger in 1926 \cite{schroedinger26} marking the beginning of wave mechanics.
+Dirac declared that chemistry has come to an end, its content being entirely contained in the powerful equation published by Schr\"odinger in 1926 \cite{schroedinger26} marking the beginning of wave mechanics.
Following the path of Schr\"odinger the problem in quantum-mechanical modeling of describing the many-body problem, i.e. a system of a large amount of interacting particles, is manifested in the high-dimensional Schr\"odinger equation for the wave function $\Psi({\vec{R}},{\vec{r}})$ that depends on the coordinates of all nuclei and electrons.
The Schr\"odinger equation contains the kinetic energy of the ions and electrons as well as the electron-ion, ion-ion and electron-electron interaction.
This cannot be solved exactly and finding approximate solutions requires several layers of simplification in order to reduce the number of free parameters.
-Approximations that consider a truncated Hilbert space of single-particle orbitals yield promising results, however, with increasing complexity and demand for high accuracy the amount of Slater determinats to be evaluated massively increases.
+Approximations that consider a truncated Hilbert space of single-particle orbitals yield promising results, however, with increasing complexity and demand for high accuracy the amount of Slater determinants to be evaluated massively increases.
In contrast, instead of using the description by the many-body wave function, the key point in density functional theory (DFT) is to recast the problem to a description utilizing the charge density $n(\vec{r})$, which constitutes a quantity in real space depending only on the three spatial coordinates.
In the following sections the basic idea of DFT will be outlined.
\text{ ,}
\end{equation}
where $Z_l$ are the atomic numbers of the nuclei and $\Psi$ is the many-electron wave function, which depends on the positions and spins of the electrons.
-Accordingly, there is only a parametrical dependence on the ionic coordinates $\vec{R}_l$.
+Accordingly, there is only a parametric dependence on the ionic coordinates $\vec{R}_l$.
However, the remaining number of free parameters is still too high and need to be further decreased.
\subsection{Hohenberg-Kohn theorem and variational principle}
This raised the question how to establish a connection between TF expressed in terms of $n(\vec{r})$ and the exact Schr\"odinger equation expressed in terms of the many-electron wave function $\Psi({\vec{r}})$ and whether a complete description in terms of the charge density is possible in principle.
The answer to this question, whether the charge density completely characterizes a system, became the starting point of modern DFT.
-Considering a system with a nondegenerate ground state there is obviously only one ground-state charge density $n_0(\vec{r})$ that correpsonds to a given potential $V(\vec{r})$.
+Considering a system with a nondegenerate ground state there is obviously only one ground-state charge density $n_0(\vec{r})$ that corresponds to a given potential $V(\vec{r})$.
In 1964 Hohenberg and Kohn showed the opposite and far less obvious result \cite{hohenberg64}.
Employing no more than the Rayleigh-Ritz minimal principle it is concluded by {\em reductio ad absurdum} that for a nondegenerate ground state the same charge density cannot be generated by different potentials.
Thus, the charge density of the ground state $n_0(\vec{r})$ uniquely determines the potential $V(\vec{r})$ and, consequently, the full Hamiltonian and ground-state energy $E_0$.
In mathematical terms the full many-electron ground state is a unique functional of the charge density.
-Im particular, $E$ is a functional $E[n(\vec{r})]$ of $n(\vec{r})$.
+In particular, $E$ is a functional $E[n(\vec{r})]$ of $n(\vec{r})$.
The ground-state charge density $n_0(\vec{r})$ minimizes the energy functional $E[n(\vec{r})]$, its value corresponding to the ground-state energy $E_0$, which enables the formulation of a minimal principle in terms of the charge density \cite{hohenberg64,levy82}
\begin{equation}
%\end{equation}
The system of interacting electrons is mapped to an auxiliary system, the Kohn-Sham (KS) system, of non-interacting electrons in an effective potential.
-The exact effective potential $V_{\text{eff}}(\vec{r})$ may be regarded as a fictious external potential yielding a gound-state density for non-interacting electrons, which is equal to that for interacting electrons in the external potential $V(\vec{r})$.
+The exact effective potential $V_{\text{eff}}(\vec{r})$ may be regarded as a fictitious external potential yielding a ground-state density for non-interacting electrons, which is equal to that for interacting electrons in the external potential $V(\vec{r})$.
The one-electron KS orbitals $\Phi_i(\vec{r})$ as well as the respective KS energies $\epsilon_i$ are not directly attached to any physical observable except for the ground-state density, which is determined by equation \eqref{eq:basics:kse3} and the ionization energy, which is equal to the highest occupied relative to the vacuum level.
The KS equations may be considered the formal exactification of the Hartree theory, which it is reduced to if the exchange-correlation potential and functional are neglected.
In addition to the Hartree-Fock (HF) method, KS theory includes the difference of the kinetic energy of interacting and non-interacting electrons as well as the remaining contributions to the correlation energy that are not part of the HF correlation.
\end{equation}
which is, thus, called local density approximation (LDA).
Here, the exchange-correlation energy per particle of the uniform electron gas of constant density $n$ is used for $\epsilon_{\text{xc}}(n(\vec{r}))$.
-Although, even in such a simple case, no exact form of the correlation part of $\epsilon_{\text{xc}}(n)$ is known, highly accurate numerical estimates using Monte Carlo methods \cite{ceperley80} and corresponding paramterizations exist \cite{perdew81}.
+Although, even in such a simple case, no exact form of the correlation part of $\epsilon_{\text{xc}}(n)$ is known, highly accurate numerical estimates using Monte Carlo methods \cite{ceperley80} and corresponding parametrizations exist \cite{perdew81}.
Obviously exact for the homogeneous electron gas, the LDA was {\em a priori} expected to be useful only for densities varying slowly on scales of the local Fermi or TF wavelength.
Nevertheless, LDA turned out to be extremely successful in describing some properties of highly inhomogeneous systems accurately within a few percent.
Although LDA is known to overestimate the cohesive energy in solids by \unit[10-20]{\%}, the lattice parameters are typically determined with an astonishing accuracy of about \unit[1]{\%}.
Finally, a set of basis functions is required to represent the one-electron KS wave functions.
With respect to the numerical treatment it is favorable to approximate the wave functions by linear combinations of a finite number of such basis functions.
-Covergence of the basis set, i.e. convergence of the wave functions with respect to the amount of basis functions, is most crucial for the accuracy of the numerical calulations.
+Convergence of the basis set, i.e. convergence of the wave functions with respect to the amount of basis functions, is most crucial for the accuracy of the numerical calculations.
Two classes of basis sets, the plane-wave and local basis sets, exist.
Local basis set functions usually are atomic orbitals, i.e. mathematical functions that describe the wave-like behavior of electrons, which are localized, i.e. centered on atoms or bonds.
The basis set is orthonormal by construction and, as mentioned above, it is simple to check for convergence.
The biggest advantage, however, is the ability to perform exact calculations by a discrete sum over a numerical grid.
This is due to the related construction of the grid and the PW basis.
-Ofcourse, exactness is restricted by the fact that the PW basis set is finite.
+Of course, exactness is restricted by the fact that the PW basis set is finite.
The simple form of the PW representation of the KS equations
\begin{equation}
\sum_{\vec{G}'} \left[
This fact is exploited in the pseudopotential (PP) approach \cite{cohen70} by removing the core electrons and replacing the atom and the associated strong ionic potential by a pseudoatom and a weaker PP that acts on a set of pseudo wave functions rather than the true valance wave functions.
This way, the pseudo wave functions become smooth near the nuclei.
-Most PPs statisfy four general conditions.
+Most PPs satisfy four general conditions.
The pseudo wave functions generated by the PP should not contain nodes, i.e. the pseudo wave functions should be smooth and free of wiggles in the core region.
Outside the core region, the pseudo and real valence wave functions as well as the generated charge densities need to be identical.
The charge enclosed within the core region must be equal for both wave functions.
V_{\text{nl}}(\vec{r}) = \sum_{lm} \mid lm \rangle V_l(\vec{r}) \langle lm \mid
\text{ .}
\end{equation}
-Applying of the operator $V_{\text{nl}}(\vec{r})$ decomposes the electronic wave functions into spherical harmonics $\mid lm \rangle$, i.e. the orbitals with azimuthal angular momentum $l$ and magnetic number $m$, which are then multiplied by the respective psuedopotential $V_l(\vec{r})$ for angular momentum $l$.
-The standard generation procedure of pseudopotentials proceeds by varying its parameters until the pseudo eigenvalues are eual to the all-electron valence eigenvalues and the pseudo wave functions match the all-electron valence wave functions beyond a certain cut-off radius detrmining the core region.
+Applying of the operator $V_{\text{nl}}(\vec{r})$ decomposes the electronic wave functions into spherical harmonics $\mid lm \rangle$, i.e. the orbitals with azimuthal angular momentum $l$ and magnetic number $m$, which are then multiplied by the respective pseudopotential $V_l(\vec{r})$ for angular momentum $l$.
+The standard generation procedure of pseudopotentials proceeds by varying its parameters until the pseudo eigenvalues are equal to the all-electron valence eigenvalues and the pseudo wave functions match the all-electron valence wave functions beyond a certain cut-off radius determining the core region.
Modified methods to generate ultra-soft pseudopotentials were proposed, which address the rapid convergence with respect to the size of the plane wave basis set \cite{vanderbilt90,troullier91}.
Using PPs the rapid oscillations of the wave functions near the core of the atoms are removed considerably reducing the number of plane waves necessary to appropriately expand the wave functions.
Since the values of the wave function within a small interval around $\vec{k}$ are almost identical, it is possible to approximate the infinite sum by a sum over an affordable number of $k$ points, each representing the respective region of the wave function in $\vec{k}$ space.
Methods have been derived for obtaining very accurate approximations by a summation over special sets of $\vec{k}$ points with distinct, associated weights \cite{baldereschi73,chadi73,monkhorst76}.
If present, symmetries in reciprocal space may further reduce the number of calculations.
-For supercells, i.e. repeating unit cells that contain several primitive cells, restricting the sampling of the Brillouin zone (BZ) to the $\Gamma$ point can yield quite accurat results.
+For supercells, i.e. repeating unit cells that contain several primitive cells, restricting the sampling of the Brillouin zone (BZ) to the $\Gamma$ point can yield quite accurate results.
In fact, with respect to BZ sampling, calculating wave functions of a supercell containing $n$ primitive cells for only one $\vec{k}$ point is equivalent to the scenario of a single primitive cell and the summation over $n$ points in $\vec{k}$ space.
-In general, finer $\vec{k}$ point meshes better account for the periodicity of a system, which in some cases, however, might be fictious anyway.
+In general, finer $\vec{k}$ point meshes better account for the periodicity of a system, which in some cases, however, might be fictitious anyway.
\subsection{Structural relaxation and Hellmann-Feynman theorem}
+\sum_j \Phi_j^* H \frac{\partial \Phi_j}{\partial \vec{R}_i}
\text{ ,}
\end{equation}
-indeed reveals a contributon to the chnage in total energy due to the change of the wave functions $\Phi_j$.
+indeed reveals a contribution to the change in total energy due to the change of the wave functions $\Phi_j$.
However, provided that the $\Phi_j$ are eigenstates of $H$, it is easy to show that the last two terms cancel each other and in the special case of $H=T+V$ the force is given by
\begin{equation}
\vec{F}_i=-\sum_j \Phi_j^*\Phi_j\frac{\partial V}{\partial \vec{R}_i}
\caption{Schematic of the constrained relaxation technique (a) and of a modified version (b) used to obtain migration pathways and corresponding configurational energies.}
\label{fig:basics:crt}
\end{figure}
-One possibility to compute the migration path from one stable cofiguration into another one is provided by the constrained relaxation technique (CRT) \cite{kaukonen98}.
+One possibility to compute the migration path from one stable configuration into another one is provided by the constrained relaxation technique (CRT) \cite{kaukonen98}.
The atoms involving great structural changes in the diffusion process are moved stepwise from the starting to the final position and relaxation after each step is only allowed in the plane perpendicular to the direction of the vector connecting its starting and final position.
This is illustrated in Fig. \ref{fig:basics:crto}.
The number of steps required for smooth transitions depends on the shape of the potential energy surface.
No constraints are applied to the remaining atoms to allow for the relaxation of the surrounding lattice.
-To prevent the remaining lattice to shift according to the displacement of the defect, ohowever, some atoms far away from the defect region should be fixed in all three coordinate directions.
+To prevent the remaining lattice to shift according to the displacement of the defect, however, some atoms far away from the defect region should be fixed in all three coordinate directions.
However, for the present study, the method tremendously failed.
-Abrupt changes in structure and configurational energy occured among relaxed structures of two successive displacement steps.
+Abrupt changes in structure and configurational energy occurred among relaxed structures of two successive displacement steps.
For some structures even the expected final configurations are not obtained.
-Thus, the method mentioned above is adjusted adding further constraints in order to obtain smooth transitions with repsect to energy and structure.
+Thus, the method mentioned above is adjusted adding further constraints in order to obtain smooth transitions with respect to energy and structure.
In the modified method all atoms are stepwise displaced towards their final positions.
In addition to this, relaxation of each individual atom is only allowed in the plane perpendicular to the last individual displacement vector, as displayed in Fig. \ref{fig:basics:crtm}.
-In the modified version respective energies could be higher than the real ones due to the additional constraints that prevent a more adequate relaxation until the final copnfiguration is reached.
+In the modified version respective energies could be higher than the real ones due to the additional constraints that prevent a more adequate relaxation until the final configuration is reached.
Structures of maximum configurational energy do not necessarily constitute saddle point configurations, i .e. the method does not guarantee to find the true minimum energy path.
-Whether a saddle point configuration and, thus, the minimum energy path is obtained by the CRT method, needs to be verified by caculating the respective vibrational modes.
+Whether a saddle point configuration and, thus, the minimum energy path is obtained by the CRT method, needs to be verified by calculating the respective vibrational modes.
Modifications used to add the CRT feature to the {\textsc vasp} code and a short instruction on how to use it can be found in appendix \ref{app:patch_vasp}.
% todo - advantages of pw basis concenring hf forces
It has turned out to be very difficult to capture the results of quantum-mechanical calculations in analytical potential models.
Among the established analytical potentials only the EDIP \cite{bazant97,justo98} and Stillinger-Weber \cite{stillinger85} potential reproduce the correct order in energy of the defects.
-However, these potenitals show shortcomings concerning the description of other physical properties and are unable to describe the C-C and C-Si interaction.
+However, these potentials show shortcomings concerning the description of other physical properties and are unable to describe the C-C and C-Si interaction.
In fact the EA potential calculations favor the tetrahedral defect configuration.
This limitation is assumed to arise due to the cut-off.
In the tetrahedral configuration the second neighbors are only slightly more distant than the first neighbors, which creates the particular problem.
While not completely rendering impossible further, more challenging empirical potential studies on large systems, the artifact has to be taken into account in the investigations of defect combinations later on in this chapter.
The hexagonal configuration is not stable opposed to results of the authors of the potential \cite{albe_sic_pot}.
-In the first two pico seconds, while kinetic energy is decoupled from the system, the \si{} seems to condense at the hexagonal site.
+In the first two picoseconds, while kinetic energy is decoupled from the system, the \si{} seems to condense at the hexagonal site.
The formation energy of \unit[4.48]{eV} is determined by this low kinetic energy configuration shortly before the relaxation process starts.
The \si{} atom then begins to slowly move towards an energetically more favorable position very close to the tetrahedral one but slightly displaced along the three coordinate axes.
The formation energy of \unit[3.96]{eV} for this type of interstitial is equal to the result for the hexagonal one in the original work \cite{albe_sic_pot}.
To exclude failures in the implementation of the potential or the MD code itself the hexagonal defect structure was double-checked with the {\textsc parcas} MD code \cite{parcas_md}.
The respective relaxation energetics are likewise plotted and look similar to the energetics obtained by {\textsc posic}.
In fact, the same type of interstitial arises using random insertions.
-In addition, variations exist, in which the displacement is only along two \hkl<1 0 0> axes ($E_\text{f}=3.8\,\text{eV}$) or along a single \hkl<1 0 0> axes ($E_\text{f}=3.6\,\text{eV}$) successively approximating the tetdrahedral configuration and formation energy.
+In addition, variations exist, in which the displacement is only along two \hkl<1 0 0> axes ($E_\text{f}=3.8\,\text{eV}$) or along a single \hkl<1 0 0> axes ($E_\text{f}=3.6\,\text{eV}$) successively approximating the tetrahedral configuration and formation energy.
The existence of these local minima located near the tetrahedral configuration seems to be an artifact of the analytical potential without physical authenticity revealing fundamental problems of analytical potential models for describing defect structures.
However, the energy barrier required for a transition into the tetrahedral configuration is small.
\begin{figure}[tp]
\end{minipage}
\end{flushleft}
\end{center}
-\caption[Relaxed C point defect configurations obtained by classical potential calculations.]{Relaxed C point defect configurations obtained by classical potential calculations. The Si/C atoms and the bonds (only for the interstitial atom) are illustrated by yellow/grey spheres and blue lines.}
+\caption[Relaxed C point defect configurations obtained by classical potential calculations.]{Relaxed C point defect configurations obtained by classical potential calculations. The Si/C atoms and the bonds (only for the interstitial atom) are illustrated by yellow/gray spheres and blue lines.}
\label{fig:defects:c_conf}
\end{figure}
\cs{} occupying an already vacant Si lattice site, which is in fact not an interstitial defect, is found to be the lowest configuration in energy for all potential models.
-An experiemntal value of the formation energy of \cs{} was determined by a fit to solubility data yielding a concentration of $3.5 \times 10^{24} \exp{(-2.3\,\text{eV}/k_{\text{B}}T)} \text{ cm}^{-3}$ \cite{bean71}.
+An experimental value of the formation energy of \cs{} was determined by a fit to solubility data yielding a concentration of $3.5 \times 10^{24} \exp{(-2.3\,\text{eV}/k_{\text{B}}T)} \text{ cm}^{-3}$ \cite{bean71}.
However, there is no particular reason for treating the prefactor as a free parameter in the fit to the experimental data.
It is simply given by the atomic density of pure silicon, which is $5\times 10^{22}\text{ cm}^{-3}$.
-Tersoff \cite{tersoff90} and Dal Pino et al. \cite{dal_pino93} pointed out that by combining this prefactor with the calculated values for the energy of formation ranging from \unit[1.6-1.89]{eV} an excellent agreement with the experimental solubility data within the entire temeprature range of the experiment is obtained.
+Tersoff \cite{tersoff90} and Dal Pino et al. \cite{dal_pino93} pointed out that by combining this prefactor with the calculated values for the energy of formation ranging from \unit[1.6-1.89]{eV} an excellent agreement with the experimental solubility data within the entire temperature range of the experiment is obtained.
This reinterpretation of the solubility data, first proposed by Tersoff and later on reinforced by Dal~Pino~et~al. is in good agreement with the results of the quantum-mechanical calculations performed in this work.
Unfortunately the EA potential undervalues the formation energy roughly by a factor of two, which is a definite drawback of the potential.
-Except for Tersoff's results for the tedrahedral configuration, the \ci{} \hkl<1 0 0> DB is the energetically most favorable interstital configuration.
+Except for Tersoff's results for the tetrahedral configuration, the \ci{} \hkl<1 0 0> DB is the energetically most favorable interstitial configuration.
As mentioned above, the low energy of formation for the tetrahedral interstitial in the case of the Tersoff potential is believed to be an artifact of the abrupt cut-off set to \unit[2.5]{\AA} (see Ref. 11 and 13 in \cite{tersoff90}) and the real formation energy is, thus, supposed to be located between \unit[3-10]{eV}.
Keeping these considerations in mind, the \ci{} \hkl<1 0 0> DB is the most favorable interstitial configuration for all interaction models.
This finding is in agreement with several theoretical\cite{burnard93,leary97,dal_pino93,capaz94,jones04} and experimental\cite{watkins76,song90} investigations, which all predict this configuration to be the ground state.
Quantum-mechanical calculations reveal this configuration to be unstable, which is also reproduced by the EA potential.
In both cases a relaxation towards the \ci{} \hkl<1 0 0> DB configuration is observed.
Opposed to results of the first-principles calculations, Tersoff finds this configuration to be stable \cite{tersoff90}.
-In fact, the stability of the hexagonal interstitial could not be reproduced in simulations performed in this work using the unmodifed Tersoff potential parameters.
+In fact, the stability of the hexagonal interstitial could not be reproduced in simulations performed in this work using the unmodified Tersoff potential parameters.
Unfortunately, apart from the modified parameters, no more conditions specifying the relaxation process are given in Tersoff's study on C point defects in Si.
The tetrahedral is the second most unfavorable interstitial configuration using classical potentials if the abrupt cut-off effect of the Tersoff potential is taken into account.
Tersoff indeed predicts a metastable BC configuration.
However, it is not in the correct order and lower in energy than the \ci{} \hkl<1 1 0> DB.
Quantum-mechanical results of this configuration are discussed in more detail in section \ref{subsection:bc}.
-In another {\em ab inito} study, Capaz~et~al.~\cite{capaz94} in turn found the BC configuration to be an intermediate saddle point structure of a possible migration path, which is \unit[2.1]{eV} higher than the \ci{} \hkl<1 0 0> DB structure.
-This is assumed to be due to the neglection of the electron spin in these calculations.
+In another {\em ab initio} study, Capaz~et~al.~\cite{capaz94} in turn found the BC configuration to be an intermediate saddle point structure of a possible migration path, which is \unit[2.1]{eV} higher than the \ci{} \hkl<1 0 0> DB structure.
+This is assumed to be due to the neglecting of the electron spin in these calculations.
Another {\textsc vasp} calculation without fully accounting for the electron spin results in the smearing of a single electron over two non-degenerate states for the BC configuration.
This problem is resolved by spin polarized calculations resulting in a net spin of one accompanied by a reduction of the total energy by \unit[0.3]{eV} and the transformation into a metastable local minimum configuration.
It is worth to note that all other listed configurations are not affected by spin polarization.
To conclude, discrepancies between the results from classical potential calculations and those obtained from first principles are observed.
Within the classical potentials EA outperforms Tersoff and is, therefore, used for further studies.
Both methods (EA and DFT) predict the \ci{} \hkl<1 0 0> DB configuration to be most stable.
-Also the remaining defects and their energetical order are described fairly well.
+Also the remaining defects and their energetic order are described fairly well.
It is thus concluded that, so far, modeling of the SiC precipitation by the EA potential might lead to trustable results.
\subsection[C \hkl<1 0 0> dumbbell interstitial configuration]{\boldmath C \hkl<1 0 0> dumbbell interstitial configuration}
\label{subsection:100db}
As the \ci{} \hkl<1 0 0> DB constitutes the ground-state configuration of a C atom incorporated into otherwise perfect c-Si it is the most probable and, hence, one of the most important interstitial configurations of C in Si.
-The structure was initially suspected by IR local vibrational mode absorption \cite{bean70} and finally verified by electron paramegnetic resonance (EPR) \cite{watkins76} studies on irradiated Si substrates at low temperatures.
+The structure was initially suspected by IR local vibrational mode absorption \cite{bean70} and finally verified by electron paramagnetic resonance (EPR) \cite{watkins76} studies on irradiated Si substrates at low temperatures.
Fig. \ref{fig:defects:100db_cmp} schematically shows the \ci{} \hkl<1 0 0> DB structure and Table \ref{tab:defects:100db_cmp} lists the details of the atomic displacements, distances and bond angles obtained by classical potential and quantum-mechanical calculations.
For comparison, the obtained structures for both methods are visualized in Fig. \ref{fig:defects:100db_vis_cmp}.
\includegraphics[height=10cm]{c_pd_vasp/eden.eps}
\includegraphics[height=12cm]{c_pd_vasp/100_2333_ksl.ps}
\end{center}
-\caption[Charge density isosurface and Kohn-Sham levels of the \ci{} \hkl<1 0 0> DB structure obtained by {\textsc vasp} calculations.]{Charge density isosurface and Kohn-Sham levels of the \ci{} \hkl<1 0 0> DB structure obtained by {\textsc vasp} calculations. Yellow and grey spheres correspond to Si and C atoms. The blue surface is the charge density isosurface. In the energy level diagram red and green lines and dots mark occupied and unoccupied states.}
+\caption[Charge density isosurface and Kohn-Sham levels of the \ci{} \hkl<1 0 0> DB structure obtained by {\textsc vasp} calculations.]{Charge density isosurface and Kohn-Sham levels of the \ci{} \hkl<1 0 0> DB structure obtained by {\textsc vasp} calculations. Yellow and gray spheres correspond to Si and C atoms. The blue surface is the charge density isosurface. In the energy level diagram red and green lines and dots mark occupied and unoccupied states.}
\label{img:defects:charge_den_and_ksl}
\end{figure}
The Si atom labeled '1' and the C atom compose the DB structure.
\caption[Structure, charge density isosurface, molecular orbital diagram and Kohn-Sham level diagram of the bond-centered interstitial configuration.]{Structure, charge density, molecular orbital diagram isosurface and Kohn-Sham level diagram of the bond-centered interstitial configuration. Gray, green and blue surfaces mark the charge density of spin up, spin down and the resulting spin up electrons in the charge density isosurface, in which the carbon atom is represented by a red sphere. In the energy level diagram red and green lines mark occupied and unoccupied states.}
\label{img:defects:bc_conf}
\end{figure}
-In the BC insterstitial configuration the interstitial atom is located in between two next neighbored Si atoms forming linear bonds.
-In a previous study this configuration was found to constitute an intermediate saddle point configuration determining the migration barrier of one possibe migration path of a \ci{} \hkl<1 0 0> DB configuration into an equivalent one \cite{capaz94}.
+In the BC interstitial configuration the interstitial atom is located in between two next neighbored Si atoms forming linear bonds.
+In a previous study this configuration was found to constitute an intermediate saddle point configuration determining the migration barrier of one possible migration path of a \ci{} \hkl<1 0 0> DB configuration into an equivalent one \cite{capaz94}.
This is in agreement with results of the EA potential simulations, which reveal this configuration to be unstable relaxing into the \ci{} \hkl<1 1 0> configuration.
However, this fact could not be reproduced by spin polarized {\textsc vasp} calculations performed in this work.
Present results suggest this configuration to correspond to a real local minimum.
-In fact, an additional barrier has to be passed to reach this configuration starting from the \ci{} \hkl<1 0 0> interstitital configuration, which is investigated in section \ref{subsection:100mig}.
+In fact, an additional barrier has to be passed to reach this configuration starting from the \ci{} \hkl<1 0 0> interstitial configuration, which is investigated in section \ref{subsection:100mig}.
After slightly displacing the C atom along the \hkl[1 0 0] (equivalent to a displacement along \hkl[0 1 0]), \hkl[0 0 1], \hkl[0 0 -1] and \hkl[1 -1 0] direction the distorted structures relax back into the BC configuration.
As will be shown in subsequent migration simulations the same would happen to structures where the C atom is displaced along the migration direction, which approximately is the \hkl[1 1 0] direction.
These relaxations indicate that the BC configuration is a real local minimum instead of an assumed saddle point configuration.
Since the activation energy of the first and last migration path is much greater than the experimental value, the second path is identified to be responsible as a migration path for the most likely C interstitial in Si explaining both, annealing and reorientation experiments.
The activation energy of roughly \unit[0.9]{eV} nicely compares to experimental values reinforcing the correct identification of the C-Si DB diffusion mechanism.
-Slightly increased values compared to experiment might be due to the tightend constraints applied in the modified CRT approach.
+Slightly increased values compared to experiment might be due to the tightened constraints applied in the modified CRT approach.
Nevertheless, the theoretical description performed in this work is improved compared to a former study \cite{capaz94}, which underestimates the experimental value by \unit[35]{\%}.
In addition, it is finally shown that the BC configuration, for which spin polarized calculations are necessary, constitutes a real local minimum instead of a saddle point configuration due to the presence of restoring forces for displacements in migration direction.
\end{figure}
Fig. \ref{fig:defects:cp_bc_00-1_mig} shows the evolution of structure and energy along the \ci{} BC to \hkl[0 0 -1] DB transition.
Since the \ci{} BC configuration is unstable relaxing into the \hkl[1 1 0] DB configuration within this potential, the low kinetic energy state is used as a starting configuration.
-Two different pathways are obtained for different time constants of the Berendse
-n thermostat.
+Two different pathways are obtained for different time constants of the Berendsen thermostat.
With a time constant of \unit[1]{fs} the C atom resides in the \hkl(1 1 0) plane
resulting in a migration barrier of \unit[2.4]{eV}.
However, weaker coupling to the heat bath realized by an increase of the time constant to \unit[100]{fs} enables the C atom to move out of the \hkl(1 1 0) plane already at the beginning, which is accompanied by a reduction in energy, approaching the final configuration on a curved path.
\end{figure}
Figures \ref{fig:defects:cp_00-1_0-10_mig} and \ref{fig:defects:cp_00-1_ip0-10_mig} show the migration barriers of the \ci{} \hkl[0 0 -1] to \hkl[0 -1 0] DB transition.
In the first case, the transition involves a change in the lattice site of the C atom whereas in the second case, a reorientation at the same lattice site takes place.
-In the first case, the pathways for the two different time cosntants look similar.
+In the first case, the pathways for the two different time constants look similar.
A local minimum exists in between two peaks of the graph.
The corresponding configuration, which is illustrated for the results obtained for a time constant of \unit[1]{fs}, looks similar to the \ci{} \hkl[1 1 0] configuration.
Indeed, this configuration is obtained by relaxation simulations without constraints of configurations near the minimum.
As mentioned earlier, the BC configuration itself constitutes a saddle point configuration relaxing into the energetically more favorable \hkl[1 1 0] DB configuration.
An activation energy of \unit[2.2]{eV} is necessary to reorientate the \hkl[0 0 -1] into the \hkl[1 1 0] DB configuration, which is \unit[1.3]{eV} higher in energy.
Residing in this state another \unit[0.90]{eV} is enough to make the C atom form a \hkl[0 0 -1] DB configuration with the Si atom of the neighbored lattice site.
-In contrast to quantum-mechanical calculations, in which the direct transition is the energetically most favorable transition and the transition composed of the intermmediate migration steps is very unlikely to occur, the just presented pathway is much more conceivable in classical potential simulations, since the energetically most favorable transition found so far is likewise composed of two migration steps with activation energies of \unit[2.2]{eV} and \unit[0.5]{eV}, for which the intermediate state is the BC configuration, which is unstable.
-Thus the just proposed migration path, which involves the \hkl[1 1 0] interstitial configuration, becomes even more probable than the initially porposed path, which involves the BC configuration that is, in fact, unstable.
+In contrast to quantum-mechanical calculations, in which the direct transition is the energetically most favorable transition and the transition composed of the intermediate migration steps is very unlikely to occur, the just presented pathway is much more conceivable in classical potential simulations, since the energetically most favorable transition found so far is likewise composed of two migration steps with activation energies of \unit[2.2]{eV} and \unit[0.5]{eV}, for which the intermediate state is the BC configuration, which is unstable.
+Thus the just proposed migration path, which involves the \hkl[1 1 0] interstitial configuration, becomes even more probable than the initially proposed path, which involves the BC configuration that is, in fact, unstable.
Due to these findings, the respective path is proposed to constitute the diffusion-describing path.
The evolution of structure and configurational energy is displayed again in Fig. \ref{fig:defects:involve110}.
\begin{figure}[tp]
The initial configuration is still evident in the relaxed configuration.
The two \ci{} atoms form a strong C-C bond, which is responsible for the large gain in energy resulting in a binding energy of \unit[-2.39]{eV}.
This bond has a length of \unit[1.38]{\AA} close to the next neighbor distance in diamond or graphite, which is approximately \unit[1.54]{\AA}.
-The minimum of the binding energy observed for this configuration suggests prefered C clustering as a competing mechnism to the \ci{} DB interstitial agglomeration inevitable for the SiC precipitation.
+The minimum of the binding energy observed for this configuration suggests preferred C clustering as a competing mechanism to the \ci{} DB interstitial agglomeration inevitable for the SiC precipitation.
However, the second most favorable configuration ($E_{\text{f}}=-2.25\,\text{eV}$) is represented four times, i.e. two times more often than the ground-state configuration, within the systematically investigated configuration space.
-Thus, particularly at high temepratures that cause an increase of the entropic contribution, this structure constitutes a serious opponent to the ground state.
+Thus, particularly at high temperatures that cause an increase of the entropic contribution, this structure constitutes a serious opponent to the ground state.
In fact, following results on migration simulations will reinforce the assumption of a low probability for C clustering by thermally activated processes.
\begin{figure}[tp]
The C atoms are spaced by \unit[3.14]{\AA}, which is very close to the expected C-C next neighbor distance of \unit[3.08]{\AA} in SiC.
Figure \ref{fig:defects:205} displays the results of a \hkl[0 0 1] DB inserted at position 3.
The binding energy is \unit[-2.05]{eV}.
-Both DBs are tilted along the same direction remaining parallely aligned and the second DB is pushed downwards in such a way, that the four DB atoms form a rhomboid.
+Both DBs are tilted along the same direction remaining aligned in parallel and the second DB is pushed downwards in such a way, that the four DB atoms form a rhomboid.
Both C atoms form tetrahedral bonds to four Si atoms.
However, Si atom number 1 and number 3, which are bound to the second \ci{} atom are also bound to the initial C atom.
These four atoms of the rhomboid reside in a plane and, thus, do not match the situation in SiC.
This low interaction is due to the large distance and a missing direct connection by bonds along a chain in the crystallographic \hkl<1 1 0> direction.
Both, C and Si atoms of the DBs form threefold coordinated bonds to their neighbors.
The energetically most unfavorable configuration ($E_{\text{b}}=0.26\,\text{eV}$) is obtained for the \ci{} \hkl[0 0 1] DB, which is oppositely orientated with respect to the initial one.
-A DB taking the same orientation as the initial one is less unfavorble ($E_{\text{b}}=0.04\,\text{eV}$).
+A DB taking the same orientation as the initial one is less unfavorable ($E_{\text{b}}=0.04\,\text{eV}$).
Both configurations are unfavorable compared to far-off, isolated DBs.
Nonparallel orientations, i.e. the \hkl[0 1 0], \hkl[0 -1 0] and its equivalents, result in binding energies of \unit[-0.12]{eV} and \unit[-0.27]{eV}, thus, constituting energetically favorable configurations.
The reduction of strain energy is higher in the second case, where the C atom of the second DB is placed in the direction pointing away from the initial C atom.
\caption[Relaxed structures of defect combinations obtained by creating {\hkl[0 0 1]}, {\hkl[0 0 -1]}, {\hkl[0 -1 0]} and {\hkl[1 0 0]} DBs at position 5.]{Relaxed structures of defect combinations obtained by creating \hkl[0 0 1] (a), \hkl[0 0 -1] (b), \hkl[0 -1 0] (c) and \hkl[1 0 0] (d) DBs at position 5.}
\label{fig:defects:comb_db_03}
\end{figure}
-Energetically beneficial configurations of defect combinations are observed for interstititals of all orientations placed at position 5, a position two bonds away from the initial interstitial along the \hkl[1 1 0] direction.
+Energetically beneficial configurations of defect combinations are observed for interstitials of all orientations placed at position 5, a position two bonds away from the initial interstitial along the \hkl[1 1 0] direction.
Relaxed structures of these combinations are displayed in Fig. \ref{fig:defects:comb_db_03}.
Fig. \ref{fig:defects:153} and \ref{fig:defects:166} show the relaxed structures of \hkl[0 0 1] and \hkl[0 0 -1] DBs.
The upper DB atoms are pushed towards each other forming fourfold coordinated bonds.
The DB structure shifts towards the position of the vacancy, which replaces the Si atom usually bound to and at the same time strained by the \si{} DB atom.
Due to the displacement into the \hkl[1 -1 0] direction the bond of the DB Si atom to the Si atom on the top left breaks and instead forms a bond to the Si atom located in \hkl[1 -1 1] direction, which is not shown in Fig.~\ref{fig:defects:314}.
A binding energy of \unit[-3.14]{eV} is obtained for this structure composing another energetically favorable configuration.
-A vacancy ctreated at position 2 enables the relaxation of Si atom number 1 mainly in \hkl[0 0 -1] direction.
+A vacancy created at position 2 enables the relaxation of Si atom number 1 mainly in \hkl[0 0 -1] direction.
The bond to Si atom number 5 breaks.
Hence, the \si{} DB atom is not only displaced along \hkl[0 0 -1] but also and to a greater extent in \hkl[1 1 0] direction.
The C atom is slightly displaced in \hkl[0 1 -1] direction.
A net amount of five Si-Si and one Si-C bond are additionally formed during transition.
An activation energy of \unit[0.6]{eV} necessary to overcome the migration barrier is found.
This energy is low enough to constitute a feasible mechanism in SiC precipitation.
-To reverse this process \unit[5.4]{eV} are needed, which make this mechanism very unprobable.
+To reverse this process \unit[5.4]{eV} are needed, which make this mechanism very improbable.
%
The migration path is best described by the reverse process.
Starting at \unit[100]{\%}, energy is needed to break the bonds of Si atom 1 to its neighbored Si atoms as well as the bond of the C atom to Si atom number 5.
With regard to the IBS process, in which highly energetic C atoms enter the Si target being able to kick out Si atoms from their lattice sites, such configurations are absolutely conceivable and a significant influence on the precipitation process might be attributed to them.
Thus, combinations of \cs{} and an additional \si{} are examined in the following.
The ground-state of a single \si{} was found to be the \si{} \hkl<1 1 0> DB configuration.
-For the follwoing study the same type of self-interstitial is assumed to provide the energetically most favorable configuration in combination with \cs.
+For the following study the same type of self-interstitial is assumed to provide the energetically most favorable configuration in combination with \cs.
\begin{table}[tp]
\begin{center}
%
In the same way the energetically most unfavorable configuration can be explained, which is configuration \RM{3}.
The \cs{} is located next to the lattice site shared by the \si{} \hkl[1 1 0] DB in \hkl[1 -1 1] direction.
-Thus, the compressive stress along \hkl[1 1 0] of the \si{} \hkl[1 1 0] DB is not compensated but intensified by the tensile stress of the \cs{} atom, which is no longer loacted along the direction of stress.
+Thus, the compressive stress along \hkl[1 1 0] of the \si{} \hkl[1 1 0] DB is not compensated but intensified by the tensile stress of the \cs{} atom, which is no longer located along the direction of stress.
However, even configuration \RM{1} is energetically less favorable than the \hkl<1 0 0> C$_{\text{i}}$ DB, which, thus, remains the ground state of a C atom introduced into otherwise perfect c-Si.
The transition involving the latter two configurations is shown in Fig.~\ref{fig:162-097}.
The possibility for separated configurations of \cs{} and \si{} becomes even more likely if one of the constituents exhibits a low barrier of migration.
In this case, the \si{} is assumed to constitute the mobile defect compared to the stable \cs{} atom.
Thus, migration paths of \si{} are investigated in the following excursus.
-Acoording to Fig.~\ref{fig:defects:si_mig1}, an activation energy of \unit[0.67]{eV} is necessary for the transition of the \si{} \hkl[0 -1 1] to \hkl[1 1 0] DB located at the neighbored Si lattice site in \hkl[1 1 -1] direction.
+According to Fig.~\ref{fig:defects:si_mig1}, an activation energy of \unit[0.67]{eV} is necessary for the transition of the \si{} \hkl[0 -1 1] to \hkl[1 1 0] DB located at the neighbored Si lattice site in \hkl[1 1 -1] direction.
\begin{figure}[tp]
\begin{center}
\includegraphics[width=0.7\textwidth]{si_110_110_mig_02_conf.ps}
However, even for non-interacting defects, the energy of formation, which is then given by the sum of the formation energies of the separated defects (\unit[4.15]{eV}) is still higher than that of the C$_{\text{i}}$ \hkl<1 0 0> DB.
Unexpectedly, the structure of a Si$_{\text{i}}$ \hkl<1 1 0> DB and a neighbored C$_{\text{s}}$, which is the most favored configuration of a C$_{\text{s}}$ and Si$_{\text{i}}$ DB according to quantum-mechanical calculations, likewise constitutes an energetically favorable configuration within the EA description, which is even preferred over the two least separated configurations of C$_{\text{s}}$ and Si$_{\text{i}}$ T.
This is attributed to an effective reduction in strain enabled by the respective combination.
-Quantum-mechanical results reveal a more favorable energy of fomation for the C$_{\text{s}}$ and Si$_{\text{i}}$ T (a) configuration.
+Quantum-mechanical results reveal a more favorable energy of formation for the C$_{\text{s}}$ and Si$_{\text{i}}$ T (a) configuration.
However, this configuration is unstable involving a structural transition into the C$_{\text{i}}$ \hkl<1 1 0> DB interstitial, thus, not maintaining the tetrahedral Si nor the \cs{} defect.
Thus, the underestimated energy of formation of C$_{\text{s}}$ within the EA calculation does not pose a serious limitation in the present context.
Nevertheless, both methods predict the same type of interstitial as the ground-state configuration and also the order in energy of the remaining defects is reproduced fairly well.
From this, a description of defect structures by classical potentials looks promising.
%
-However, focussing on the description of diffusion processes the situation changes completely.
+However, focusing on the description of diffusion processes the situation changes completely.
Qualitative and quantitative differences exist.
First of all, a different pathway is suggested as the lowest energy path, which again might be attributed to the absence of quantum-mechanical effects in the classical interaction model.
Secondly, the activation energy is overestimated by a factor of 2.4 to 3.5 compared to the more accurate quantum-mechanical methods and experimental findings.
Quantum-mechanical investigations of two \ci{} of the \hkl<1 0 0>-type and varying separations and orientations state an attractive interaction between these interstitials.
Obtained results for the most part compare well with results gained in previous studies \cite{leary97,capaz98,mattoni2002,liu02} and show an astonishingly good agreement with experiment \cite{song90}.
%
-Depending on orientation, energetically favorable configurations are found, in which these two interstitials are located close together instead of the occurernce of largely separated and isolated defects.
+Depending on orientation, energetically favorable configurations are found, in which these two interstitials are located close together instead of the occurrence of largely separated and isolated defects.
This is due to strain compensation enabled by the combination of such defects in certain orientations.
For dumbbells oriented along the \hkl<1 1 0> bond chain and the assumption that there is the possibility of free orientation, an interaction energy proportional to the reciprocal cube of the distance in the far field regime is found.
These findings support the assumption of the \ci{} DB agglomeration.
Investigated combinations of C interstitials, however, result in distorted configurations, in which C atoms, which at some point will form SiC, are no longer aligned to the host.
-It is easily understandable that the mismatch in alignement will increase with increasing defect density.
+It is easily understandable that the mismatch in alignment will increase with increasing defect density.
In addition, the amount of Si self-interstitials equal to the amount of agglomerated C atoms would be released all of a sudden probably not being able to diffuse into the c-Si host matrix without damaging the Si surrounding or the precipitate itself.
In addition, IBS results in the formation of the cubic polytype of SiC only.
-As this result conforms well with the model of precipitation by substitutional C there is no obvious reason why hexagonal polytypes should not be able to form or an equal alignement would be mandatory assuming the model of precipitation by C-Si dumbbell agglomeration.
+As this result conforms well with the model of precipitation by substitutional C there is no obvious reason why hexagonal polytypes should not be able to form or an equal alignment would be mandatory assuming the model of precipitation by C-Si dumbbell agglomeration.
\fi
Silicon carbide (SiC) has a number of remarkable physical and chemical properties that make it a promising new material in various fields of applications.
The high electron mobility and saturation drift velocity as well as the high band gap and breakdown field in conjunction with its unique thermal stability and conductivity unveil SiC as the ideal candidate for high-power, high-frequency and high-temperature electronic and optoelectronic devices exceeding conventional silicon based solutions \cite{wesch96,morkoc94,casady96,capano97,pensl93}.
-Due to the large Si--C bonding energy SiC is a hard and chemical inert material suitable for applications under extreme conditions and capable for microelectromechanical systemis, both as structural material and as a coating layer \cite{sarro00,park98}.
+Due to the large Si--C bonding energy SiC is a hard and chemical inert material suitable for applications under extreme conditions and capable for microelectromechanical systems, both as structural material and as a coating layer \cite{sarro00,park98}.
Its radiation hardness allows the operation as a first wall material in nuclear reactors \cite{giancarli98} and as electronic devices in space \cite{capano97}.
The realization of silicon carbide based applications demands for reasonable sized wafers of high crystalline quality.
-Despite the tremendous progress achieved in the fabrication of high purity SiC employing techniques like the modified Lely process for bulk crystal growth \cite{tairov78,tsvetkov98} or chemical vapour deposition (CVD) and molecular beam epitaxy (MBE) for homo- and heteroepitaxial growth \cite{kimoto93,powell90,fissel95}, available wafer dimensions and crystal qualities are not yet considered sufficient enough.
+Despite the tremendous progress achieved in the fabrication of high purity SiC employing techniques like the modified Lely process for bulk crystal growth \cite{tairov78,tsvetkov98} or chemical vapor deposition (CVD) and molecular beam epitaxy (MBE) for homo- and heteroepitaxial growth \cite{kimoto93,powell90,fissel95}, available wafer dimensions and crystal qualities are not yet considered sufficient enough.
Another promising alternative to fabricate SiC is ion beam synthesis (IBS).
High-dose carbon implantation at elevated temperatures into silicon with subsequent annealing results in the formation of buried epitaxial SiC layers \cite{borders71,reeson87}.
-A two-temperature implantation technique was proposed to achieve single crytalline SiC layers and a sharp SiC/Si interface \cite{lindner99,lindner99_2,lindner01,lindner02}.
+A two-temperature implantation technique was proposed to achieve single crystalline SiC layers and a sharp SiC/Si interface \cite{lindner99,lindner99_2,lindner01,lindner02}.
Although high-quality SiC can be achieved by means of IBS the precipitation mechanism is not yet fully understood.
-High resolution transmisson electron microscopy studies indicate the formation of C-Si interstitial complexes sharing conventional silicon lattice sites (C-Si dumbbells) during the implantation of carbon in silicon.
+High resolution transmission electron microscopy studies indicate the formation of C-Si interstitial complexes sharing conventional silicon lattice sites (C-Si dumbbells) during the implantation of carbon in silicon.
These C-Si dumbbells agglomerate and once a critical radius is reached, the topotactic transformation into a SiC precipitate occurs \cite{werner97,lindner01}.
In contrast, investigations of strained Si$_{1-y}$C$_y$/Si heterostructures form
ed by MBE~\cite{strane94,guedj98}, which incidentally involve the formation of SiC nanocrystallites, suggest an initial coherent precipitation by agglomeration of substitutional instead of interstitial C.
In chapter \ref{chapter:sic_rev} a review of the Si/C compound is given including the very central discussion on two controversial precipitation mechanisms present in literature in section \ref{section:assumed_prec}.
Chapter \ref{chapter:basics} introduces some basics and internals of the utilized atomistic simulations as well as special methods of application.
Details of the simulation and associated test calculations are presented in chapter \ref{chapter:simulation}.
-In chapter \ref{chapter:defects} results of investigations of single defect configurations, structures of comnbinations of two individual defects as well as some selected diffusion pathways in silicon are shown.
+In chapter \ref{chapter:defects} results of investigations of single defect configurations, structures of combinations of two individual defects as well as some selected diffusion pathways in silicon are shown.
These allow to draw conclusions with respect to the SiC precipitation mechanism in Si.
More complex systems aiming to model the transformation of C incorporated in bulk Si into a SiC nucleus are examined in chapter \ref{chapter:md}.
Finally a summary and some concluding remarks are given in chapter \ref{chapter:summary}.
At the beginning, simulations are performed, which try to mimic the conditions during IBS.
Results reveal limitations of the employed potential and MD in general.
With reference to the results of the last chapter, a workaround is discussed.
-The approach is follwed and, finally, results gained by the MD simulations are interpreted drawing special attention to the established controversy concerning precipitation of SiC in Si.
+The approach is followed and, finally, results gained by the MD simulations are interpreted drawing special attention to the established controversy concerning precipitation of SiC in Si.
\section{Simulations at temperatures used in IBS}
\label{section:initial_sims}
\hline
\end{tabular}
\end{center}
-\caption{Side lengthes of the insertion volumes $V_1$, $V_2$ and $V_3$ used for the incoorperation of 6000 C atoms.}
+\caption{Side lengths of the insertion volumes $V_1$, $V_2$ and $V_3$ used for the incorporation of 6000 C atoms.}
\label{table:md:ins_vols}
\end{table}
\begin{pspicture}(0,0)(15,17)
\psframe*[linecolor=hb](3,11.5)(11,17)
- \rput[lt](3.2,16.8){\color{gray}INITIALIZIATION}
+ \rput[lt](3.2,16.8){\color{gray}INITIALIZATION}
\rput(7,16){\rnode{14}{\psframebox{Create $31\times 31\times 31$
unit cells of c-Si}}}
\rput(7,15){\rnode{13}{\psframebox{$T_{\text{s}}=450\,^{\circ}\mathrm{C}$,
\psframe*[linecolor=lbb](3,6.5)(11,11)
\rput[lt](3.2,10.8){\color{gray}CARBON INSERTION}
\rput(3,10.8){\pnode{CI}}
- \rput(7,10){\rnode{9}{\psframebox{Insertion of 10 C aoms}}}
+ \rput(7,10){\rnode{9}{\psframebox{Insertion of 10 C atoms}}}
\rput(7,9){\rnode{8}{\psframebox{Continue for 100 fs}}}
\rput(7,8){\rnode{7}{\psframebox{$T_{\text{avg}}=T_{\text{s}}
\pm1\,^{\circ}\mathrm{C}$}}}
\psset{fillcolor=white}
\nbput*{\scriptsize true}
- \rput(7,5.75){\rnode{5}{\psframebox{Continue for 100 ps}}}
+ \rput(7,5.75){\rnode{5}{\psframebox{Continue for \unit[100]{ps}}}}
\ncline[]{->}{6}{5}
\trput*{\scriptsize true}
\end{figure}
It is easily and instantly visible that there is no significant difference among the two simulations of high C concentration.
Thus, in the following, the focus can indeed be directed to low ($V_1$) and high ($V_2$, $V_3$) C concentration simulations.
-The first C-C peak appears at about \unit[0.15]{nm}, which is compareable to the nearest neighbor distance of graphite or diamond.
+The first C-C peak appears at about \unit[0.15]{nm}, which is comparable to the nearest neighbor distance of graphite or diamond.
The number of C-C bonds is much smaller for $V_1$ than for $V_2$ and $V_3$ since C atoms are spread over the total simulation volume.
On average, there are only 0.2 C atoms per Si unit cell.
These C atoms are assumed to form strong bonds.
However, no clear peak is observed but the interval of enhanced $g(r)$ values corresponds to the width of the C-C $g(r)$ peak.
In addition the abrupt increase of Si pairs at \unit[0.29]{nm} can be attributed to the Si-Si cut-off radius of \unit[0.296]{nm} as used in the present bond order potential.
The cut-off function causes artificial forces pushing the Si atoms out of the cut-off region.
-Without the abrubt increase, a maximum around \unit[0.31]{nm} gets even more conceivable.
-Analyses of randomly chosen configurations, in which distances around \unit[0.3]{nm} appear, identify \ci{} \hkl<1 0 0> DBs to be responsible for stretching the Si-Si next neighbour distance for low C concentrations, i.e. for the $V_1$ and early stages of $V_2$ and $V_3$ simulation runs.
+Without the abrupt increase, a maximum around \unit[0.31]{nm} gets even more conceivable.
+Analyses of randomly chosen configurations, in which distances around \unit[0.3]{nm} appear, identify \ci{} \hkl<1 0 0> DBs to be responsible for stretching the Si-Si next neighbor distance for low C concentrations, i.e. for the $V_1$ and early stages of $V_2$ and $V_3$ simulation runs.
This excellently agrees with the calculated value $r(13)$ in Table~\ref{tab:defects:100db_cmp} for a resulting Si-Si distance in the \ci \hkl<1 0 0> DB configuration.
\begin{figure}[tp]
The first peak observed for all insertion volumes is at approximately \unit[0.186]{nm}.
This corresponds quite well to the expected next neighbor distance of \unit[0.189]{nm} for Si and C atoms in 3C-SiC.
By comparing the resulting Si-C bonds of a \ci{} \hkl<1 0 0> DB with the C-Si distances of the low concentration simulation, it is evident that the resulting structure of the $V_1$ simulation is clearly dominated by this type of defect.
-This is not surpsising, since the \ci{} \hkl<1 0 0> DB is found to be the ground-state defect of a C interstitial in c-Si and, for the low concentration simulations, a C interstitial is expected in every fifth Si unit cell only, thus, excluding defect superposition phenomena.
+This is not surprising, since the \ci{} \hkl<1 0 0> DB is found to be the ground-state defect of a C interstitial in c-Si and, for the low concentration simulations, a C interstitial is expected in every fifth Si unit cell only, thus, excluding defect superposition phenomena.
The peak distance at \unit[0.186]{nm} and the bump at \unit[0.175]{nm} corresponds to the distance $r(3C)$ and $r(1C)$ as listed in Table~\ref{tab:defects:100db_cmp} and visualized in Fig.~\ref{fig:defects:100db_cmp}.
In addition, it can be easily identified that the \ci{} \hkl<1 0 0> DB configuration contributes to the peaks at about \unit[0.335]{nm}, \unit[0.386]{nm}, \unit[0.434]{nm}, \unit[0.469]{nm} and \unit[0.546]{nm} observed in the $V_1$ simulation.
Not only the peak locations but also the peak widths and heights become comprehensible.
The distinct peak at \unit[0.26]{nm}, which exactly matches the cut-off radius of the Si-C interaction, is again a potential artifact.
-For high C concentrations, i.e. the $V_2$ and $V_3$ simulation corresponding to a C density of about 8 atoms per c-Si unit cell, the defect concentration is likewiese increased and a considerable amount of damage is introduced in the insertion volume.
+For high C concentrations, i.e. the $V_2$ and $V_3$ simulation corresponding to a C density of about 8 atoms per c-Si unit cell, the defect concentration is likewise increased and a considerable amount of damage is introduced in the insertion volume.
The consequential superposition of these defects and the high amounts of damage generate new displacement arrangements for the C-C as well as for the Si-C pair distances, which become hard to categorize and trace and obviously lead to a broader distribution.
Short range order indeed is observed, i.e. the large amount of strong neighbored C-C bonds at \unit[0.15]{nm} as expected in graphite or diamond and Si-C bonds at \unit[0.19]{nm} as expected in SiC, but only hardly visible is the long range order.
This indicates the formation of an amorphous SiC-like phase.
The first reason is a general problem of MD simulations in conjunction with limitations in computer power, which results in a slow and restricted propagation in phase space.
In molecular systems, characteristic motions take place over a wide range of time scales.
Vibrations of the covalent bond take place on the order of \unit[10$^{-14}$]{s}, of which the thermodynamic and kinetic properties are well described by MD simulations.
-To avoid dicretization errors, the integration timestep needs to be chosen smaller than the fastest vibrational frequency in the system.
+To avoid discretization errors, the integration time step needs to be chosen smaller than the fastest vibrational frequency in the system.
On the other hand, infrequent processes, such as conformational changes, reorganization processes during film growth, defect diffusion and phase transitions are processes undergoing long-term evolution in the range of microseconds.
This is due to the existence of several local minima in the free energy surface separated by large energy barriers compared to the kinetic energy of the particles, i.e. the system temperature.
-Thus, the average time of a transition from one potential basin to another corresponds to a great deal of vibrational periods, which in turn determine the integration timestep.
+Thus, the average time of a transition from one potential basin to another corresponds to a great deal of vibrational periods, which in turn determine the integration time step.
Hence, time scales covering the necessary amount of infrequent events to observe long-term evolution are not accessible by traditional MD simulations, which are limited to the order of nanoseconds.
New methods have been developed to bypass the time scale problem.
-The most famous appraiches are hyperdnyamics (HMD) \cite{voter97,voter97_2}, parallel replica dynamics \cite{voter98}, temperature acclerated dynamics (TAD) \cite{sorensen2000} and self-guided dynamics (SGMD) \cite{wu99}, which accelerate phase space propagation while retaining proper thermodynmic sampling.
+The most famous approaches are hyperdynamics (HMD) \cite{voter97,voter97_2}, parallel replica dynamics \cite{voter98}, temperature accelerated dynamics (TAD) \cite{sorensen2000} and self-guided dynamics (SGMD) \cite{wu99}, which accelerate phase space propagation while retaining proper thermodynamic sampling.
In addition to the time scale limitation, problems attributed to the short range potential exist.
-The sharp cut-off funtion, which limits the interacting ions to the next neighbored atoms by gradually pushing the interaction force and energy to zero between the first and second next neighbor distance, is responsible for overestimated and unphysical high forces of next neighbored atoms \cite{tang95,mattoni2007}.
+The sharp cut-off function, which limits the interacting ions to the next neighbored atoms by gradually pushing the interaction force and energy to zero between the first and second next neighbor distance, is responsible for overestimated and unphysical high forces of next neighbored atoms \cite{tang95,mattoni2007}.
This is supported by the overestimated activation energies necessary for C diffusion as investigated in section \ref{subsection:defects:mig_classical}.
-Indeed, it is not only the strong C-C bond, which is hard to break, inhibiting C diffusion and further rearrengements.
+Indeed, it is not only the strong C-C bond, which is hard to break, inhibiting C diffusion and further rearrangements.
This is also true for the low concentration simulations dominated by the occurrence of C-Si DBs spread over the whole simulation volume.
The bonds of these C-Si pairs are also affected by the cut-off artifact preventing C diffusion and agglomeration of the DBs.
This can be seen from the almost horizontal progress of the total energy graph in the continuation step of Fig.~\ref{fig:md:energy_450}, even for the low concentration simulation.
These unphysical effects inherent to this type of model potentials are solely attributed to their short range character.
While cohesive and formational energies are often well described, these effects increase for non-equilibrium structures and dynamics.
-However, since valueable insights into various physical properties can be gained using this potentials, modifications mainly affecting the cut-off were designed.
+However, since valuable insights into various physical properties can be gained using this potentials, modifications mainly affecting the cut-off were designed.
One possibility is to simply skip the force contributions containing the derivatives of the cut-off function, which was successfully applied to reproduce the brittle propagation of fracture in SiC at zero temperature \cite{mattoni2007}.
-Another one is to use variable cut-off values scaled by the system volume, which properly describes thermomechanical properties of 3C-SiC \cite{tang95} but might be rather ineffective for the challange inherent to this study.
+Another one is to use variable cut-off values scaled by the system volume, which properly describes thermomechanical properties of 3C-SiC \cite{tang95} but might be rather ineffective for the challenge inherent to this study.
To conclude the obstacle needed to get passed is twofold.
The sharp cut-off of the employed bond order model potential introduces overestimated high forces between next neighbored atoms enhancing the problem of slow phase space propagation immanent to MD simulations.
This obstacle could be referred to as {\em potential enhanced slow phase space propagation}.
-Due to this, pushing the time scale to the limits of computational ressources or applying one of the above mentioned accelerated dynamics methods exclusively will not be sufficient enough.
+Due to this, pushing the time scale to the limits of computational resources or applying one of the above mentioned accelerated dynamics methods exclusively will not be sufficient enough.
Instead, the approach followed in this study, is the use of higher temperatures as exploited in TAD to find transition pathways of one local energy minimum to another one more quickly.
Since merely increasing the temperature leads to different equilibrium kinetics than valid at low temperatures, TAD introduces basin-constrained MD allowing only those transitions that should occur at the original temperature and a properly advancing system clock \cite{sorensen2000}.
This is justified by two reasons.
First of all, a compensation of the overestimated bond strengths due to the short range potential is expected.
Secondly, there is no conflict applying higher temperatures without the TAD corrections, since crystalline 3C-SiC is also observed for higher temperatures than \unit[450]{$^{\circ}$C} in IBS \cite{nejim95,lindner01}.
-It is therefore expected that the kinetics affecting the 3C-SiC precipitation are not much different at higher temperatures aside from the fact that it is occuring much more faster.
+It is therefore expected that the kinetics affecting the 3C-SiC precipitation are not much different at higher temperatures aside from the fact that it is occurring much more faster.
Moreover, the interest of this study is focused on structural evolution of a system far from equilibrium instead of equilibrium properties which rely upon proper phase space sampling.
-On the other hand, during implantation, the actual temperature inside the implantation volume is definetly higher than the experimentally determined temperature tapped from the surface of the sample.
+On the other hand, during implantation, the actual temperature inside the implantation volume is definitely higher than the experimentally determined temperature tapped from the surface of the sample.
\section{Increased temperature simulations}
\label{section:md:inct}
\label{eq:md:qdef}
\end{equation}
By this, bulk 3C-SiC will still result in $Q=1$ and precipitates will also reach values close to one.
-However, since the quality value does not account for bond lengthes, bond angles, crystallinity or the stacking sequence, high values of $Q$ not necessarily correspond to structures close to 3C-SiC.
+However, since the quality value does not account for bond lengths, bond angles, crystallinity or the stacking sequence, high values of $Q$ not necessarily correspond to structures close to 3C-SiC.
Structures that look promising due to high quality values need to be further investigated by other means.
\subsection{Low C concentration simulations}
\includegraphics[width=0.7\textwidth]{tot_pc_thesis.ps}\\
\includegraphics[width=0.7\textwidth]{tot_ba.ps}
\end{center}
-\caption[Si-C radial distribution and evolution of quality $Q$ for the low concentration simulations at different elevated temperatures.]{Si-C radial distribution and evolution of quality $Q$ according to equation \ref{eq:md:qdef} for the low concentration simulations at different elevated temperatures. All structures are cooled down to \degc{20}. The grey line shows resulting Si-C bonds in a configuration of \cs{} in c-Si (C$_\text{sub}$) at zero temperature. Arrows in the quality plot mark the end of C insertion and the start of the cooling down step. A fit function according to equation \eqref{eq:md:fit} shows the estimated evolution of quality in the absence of the cooling down sequence.}
+\caption[Si-C radial distribution and evolution of quality $Q$ for the low concentration simulations at different elevated temperatures.]{Si-C radial distribution and evolution of quality $Q$ according to equation \ref{eq:md:qdef} for the low concentration simulations at different elevated temperatures. All structures are cooled down to \degc{20}. The gray line shows resulting Si-C bonds in a configuration of \cs{} in c-Si (C$_\text{sub}$) at zero temperature. Arrows in the quality plot mark the end of C insertion and the start of the cooling down step. A fit function according to equation \eqref{eq:md:fit} shows the estimated evolution of quality in the absence of the cooling down sequence.}
\label{fig:md:tot_si-c_q}
\end{figure}
-Fig.~\ref{fig:md:tot_si-c_q} shows the radial distribution of Si-C bonds for different temperatures and the corresponding evolution of quality $Q$ as defined above for the low concentration simulaton.
+Fig.~\ref{fig:md:tot_si-c_q} shows the radial distribution of Si-C bonds for different temperatures and the corresponding evolution of quality $Q$ as defined above for the low concentration simulation.
The first noticeable and promising change in the Si-C radial distribution is the successive decline of the artificial peak at the Si-C cut-off distance with increasing temperature up to the point of disappearance at temperatures above \degc{1650}.
Obviously, sufficient kinetic energy is provided to affected atoms that are enabled to escape the cut-off region.
Additionally, a more important structural change is observed, which is illustrated in the two shaded areas in Fig.~\ref{fig:md:tot_si-c_q}.
%
-In the grey shaded region a decrease of the peak at \unit[0.186]{nm} and of the bump at \distn{0.175} accompanied by an increase of the peak at \distn{0.197} with increasing temperature is visible.
+In the gray shaded region a decrease of the peak at \unit[0.186]{nm} and of the bump at \distn{0.175} accompanied by an increase of the peak at \distn{0.197} with increasing temperature is visible.
Similarly, the peaks at \distn{0.335} and \distn{0.386} shrink in contrast to a new peak forming at \distn{0.372} as can be seen in the yellow shaded region.
Obviously, the structure obtained from the \degc{450} simulations, which is dominated by the existence of \ci{} \hkl<1 0 0> DBs, transforms into a different structure with increasing simulation temperature.
Investigations of the atomic data reveal \cs{} to be responsible for the new Si-C bonds.
While simulations at \degc{450} exhibit \perc{10} of fourfold coordinated C, simulations at \degc{2050} exceed the \perc{80} range.
Since \cs{} has four nearest neighbored Si atoms and is the preferential type of defect in elevated temperature simulations, the increase of the quality values become evident.
The quality values at a fixed temperature increase with simulation time.
-After the end of the insertion sequence marked by the first arrow, the quality is increasing and a saturation behaviour, yet before the cooling process starts, can be expected.
+After the end of the insertion sequence marked by the first arrow, the quality is increasing and a saturation behavior, yet before the cooling process starts, can be expected.
The evolution of the quality value of the simulation at \degc{2050} inside the range, in which the simulation is continued at constant temperature for \unit[100]{fs}, is well approximated by the simple fit function
\begin{equation}
f(t)=a-\frac{b}{t} \text{ ,}
\label{fig:md:tot_c-c_si-si}
\end{figure}
The formation of \cs{} also affects the Si-Si radial distribution displayed in the lower part of Fig.~\ref{fig:md:tot_c-c_si-si}.
-Investigating the atomic strcuture indeed shows that the peak arising at \distn{0.325} with increasing temperature is due to two Si atoms that form direct bonds to the \cs{} atom.
+Investigating the atomic structure indeed shows that the peak arising at \distn{0.325} with increasing temperature is due to two Si atoms that form direct bonds to the \cs{} atom.
The peak corresponds to the distance of next neighbored Si atoms along the \hkl<1 1 0> bond chain with C$_{\text{s}}$ in between.
Since the expected distance of these Si pairs in 3C-SiC is \distn{0.308}, the existing SiC structures embedded in the c-Si host are stretched.
In the upper part of Fig.~\ref{fig:md:tot_c-c_si-si} the C-C radial distribution is shown.
The total amount of C-C bonds with a distance inside the displayed separation range does not change significantly.
-Thus, even for elevated temperatures, agglomeration of C atoms neccessary to form a SiC precipitate does not take place within the simulated time scale.
+Thus, even for elevated temperatures, agglomeration of C atoms necessary to form a SiC precipitate does not take place within the simulated time scale.
However, with increasing temperature, a decrease of the amount of next neighbored C pairs can be observed.
-This is a promising result gained by the high-temperature simulations since the breaking of these diomand and graphite like bonds is mandatory for the formation of 3C-SiC.
+This is a promising result gained by the high-temperature simulations since the breaking of these diamond and graphite like bonds is mandatory for the formation of 3C-SiC.
Obviously, acceleration of the dynamics occurred by supplying additional kinetic energy.
A slight shift towards higher distances can be observed for the maximum located shortly above \distn{0.3}.
Arrows with dashed lines mark C-C distances resulting from \ci{} \hkl<1 0 0> DB combinations while arrows with solid lines mark distances arising from combinations of \cs.
%However, this is contrary to the initial precipitation model proposed in section \ref{section:assumed_prec}, which assumes that the transformation into 3C-SiC takes place in a very last step once enough C-Si DBs agglomerated.
To summarize, results of low concentration simulations show a phase transition in conjunction with an increase in temperature.
-The \ci{} \hkl<1 0 0> DB dominated struture turns into a structure characterized by the occurence of an increasing amount of \cs{} with respect to temperature.
+The \ci{} \hkl<1 0 0> DB dominated structure turns into a structure characterized by the occurrence of an increasing amount of \cs{} with respect to temperature.
Clearly, the high-temperature results indicate the precipitation mechanism involving an increased participation of \cs.
Although diamond and graphite like bonds are reduced, no agglomeration of C is observed within the simulated time.
-Isolated structures of stretched SiC, which are adjusted to the c-Si host with respect to the lattice constant and alignement, are formed.
+Isolated structures of stretched SiC, which are adjusted to the c-Si host with respect to the lattice constant and alignment, are formed.
By agglomeration of \cs{} the interfacial energy could be overcome and a transition from a coherent and stretched SiC structure into an incoherent and partially strain-compensated SiC precipitate could occur.
Indeed, \si in the near surrounding is observed, which may initially compensate tensile strain in the stretched SiC structure or rearrange the \cs{} sublattice and finally serve as supply for additional C to form further SiC or compensate strain at the interface of the incoherent SiC precipitate and the Si host.
Again, in both cases, the cut-off artifact decreases with increasing temperature.
Peaks that already exist for the low temperature simulations get slightly more distinct for elevated temperatures.
This is also true for peaks located past distances of next neighbors indicating an increase in the long range order.
-However, this change is rather small and no significant structural change is observeable.
+However, this change is rather small and no significant structural change is observable.
Due to the continuity of high amounts of damage, atomic configurations remain hard to identify even for the highest temperature.
-Other than in the low concentration simulation, analyzed defect structures are no longer necessarily aligned to the primarily existing but succesively disappearing c-Si host matrix inhibiting or at least hampering their identification and classification.
+Other than in the low concentration simulation, analyzed defect structures are no longer necessarily aligned to the primarily existing but successively disappearing c-Si host matrix inhibiting or at least hampering their identification and classification.
As for low temperatures, order in the short range exists decreasing with increasing separation.
The increase of the amount of Si-C pairs at \distn{0.186} could be positively interpreted since this type of bond also exists in 3C-SiC.
On the other hand, the amount of next neighbored C atoms with a distance of approximately \distn{0.15}, which is the distance of C in graphite or diamond, is likewise increased.
Thus, higher temperatures seem to additionally enhance a conflictive process, i.e. the formation of C agglomerates, obviously inconsistent with the desired process of 3C-SiC formation.
This is supported by the C-C peak at \distn{0.252}, which corresponds to the second next neighbor distance in the diamond structure of elemental C.
Investigating the atomic data indeed reveals two C atoms, which are bound to and interconnected by a third C atom, to be responsible for this distance.
-The C-C peak at about \distn{0.31}, wich is slightly shifted to higher distances (\distn{0.317}) with increasing temperature still corresponds quite well to the next neighbor distance of C in 3C-SiC as well as a-SiC and indeed results from C-Si-C bonds.
+The C-C peak at about \distn{0.31}, which is slightly shifted to higher distances (\distn{0.317}) with increasing temperature still corresponds quite well to the next neighbor distance of C in 3C-SiC as well as a-SiC and indeed results from C-Si-C bonds.
The Si-C peak at \distn{0.282}, which is pronounced with increasing temperature, is constructed out of a Si atom and a C atom, which are both bound to another central C atom.
This is similar for the Si-C peak at approximately \distn{0.35}.
In this case, the Si and the C atom are bound to a central Si atom.
Thus, for instance, it is not surprising that short range potentials show overestimated melting temperatures while properties of structures that are only slightly deviated from equilibrium are well described.
Due to this, increased temperatures are considered exceptionally necessary for modeling non-equilibrium processes and structures such as IBS and 3C-SiC.
-Thus, it is concluded that increased temperatures are not exclusively usefull to accelerate the dynamics approximatively describing the structural evolution.
+Thus, it is concluded that increased temperatures are not exclusively useful to accelerate the dynamics approximatively describing the structural evolution.
Moreover, it can be considered a necessary condition to deviate the system out of equilibrium enabling the formation of 3C-SiC, which is obviously realized by a successive agglomeration of \cs{}.
\ifnum1=0
As discussed in section~\ref{section:md:limit} and~\ref{section:md:inct} a further increase of the system temperature might help to overcome limitations of the short range potential and accelerate the dynamics involved in structural evolution.
Furthermore, these results indicate that increased temperatures are necessary to drive the system out of equilibrium enabling conditions needed for the formation of a metastable cubic polytype of SiC.
-A maximum temperature to avoid melting is determined in section \ref{section:md:tval} to be 120 \% of the Si melting point but due to defects lowering the transition point a maximum temperature of 95 \% of the Si melting temperature is considered usefull.
+A maximum temperature to avoid melting is determined in section \ref{section:md:tval} to be 120 \% of the Si melting point but due to defects lowering the transition point a maximum temperature of 95 \% of the Si melting temperature is considered useful.
This value is almost equal to the temperature of $2050\,^{\circ}\mathrm{C}$ already used in former simulations.
Since the maximum temperature is reached the approach is reduced to the application of longer time scales.
-This is considered usefull since the estimated evolution of quality in the absence of the cooling down sequence in figure~\ref{fig:md:tot_si-c_q} predicts an increase in quality and, thus, structural evolution is liekyl to occur if the simulation is proceeded at maximum temperature.
+This is considered useful since the estimated evolution of quality in the absence of the cooling down sequence in figure~\ref{fig:md:tot_si-c_q} predicts an increase in quality and, thus, structural evolution is likely to occur if the simulation is proceeded at maximum temperature.
Next to the employment of longer time scales and a maximum temperature a few more changes are applied.
In the following simulations the system volume, the amount of C atoms inserted and the shape of the insertion volume are modified from the values used in first MD simulations.
To speed up the simulation the initial simulation volume is reduced to 21 Si unit cells in each direction and 5500 inserted C atoms in either the whole volume or in a sphere with a radius of 3 nm corresponding to the size of a precipitate consisting of 5500 C atoms.
-The 100 ps sequence after C insertion intended for structural evolution is exchanged by a 10 ns sequence, which is hoped to result in the occurence of infrequent processes and a subsequent phase transition.
-The return to lower temperatures is considered seperately.
+The \unit[100]{ps} sequence after C insertion intended for structural evolution is exchanged by a \unit[10]{ns} sequence, which is hoped to result in the occurrence of infrequent processes and a subsequent phase transition.
+The return to lower temperatures is considered separately.
\begin{figure}[tp]
\begin{center}
\end{figure}
Fig.~\ref{fig:md:95_long_time_v1} shows the evolution in time of the radial distribution for Si-C and C-C pairs for a low C concentration simulation.
-Differences are observed for both types of atom pairs indeed indicating proceeding structural changes even well beyond 100 ps of simulation time.
+Differences are observed for both types of atom pairs indeed indicating proceeding structural changes even well beyond \unit[100]{ps} of simulation time.
Peaks attributed to the existence of substitutional C increase and become more distinct.
This finding complies with the predicted increase of quality evolution as explained earlier.
More and more C forms tetrahedral bonds to four Si neighbors occupying vacant Si sites.
There are only small changes identifiable.
A slight increase of the Si-C peak at approximately 0.36 nm attributed to the distance of substitutional C and the next but one Si atom along \hkl<1 1 0> is observed.
In the same time the C-C peak at approximately 0.32 nm corresponding to the distance of two C atoms interconnected by a Si atom along \hkl<1 1 0> slightly decreases.
-Obviously the system preferes a slight increase of isolated substitutional C at the expense of incoherent C-Si-C precipitate configurations, which at a first glance actually appear as promising configurations in the precipitation event.
+Obviously, the system prefers a slight increase of isolated substitutional C at the expense of incoherent C-Si-C precipitate configurations, which at a first glance actually appear as promising configurations in the precipitation event.
On second thoughts however, this process of splitting a C atom out of this structure is considered necessary in order to allow for the rearrangement of C atoms on substitutional lattice sites on the one hand and for C diffusion otherwise, which is needed to end up in a structure, in which one of the two fcc sublattices is composed out of C only.
For both, high and low concentration simulations the radial distribution converges as can be seen by the nearly identical graphs of the two most advanced configurations.
For the low C concentrations, time scales are still too low to observe C agglomeration sufficient for SiC precipitation, which is attributed to the slow phase space propagation inherent to MD in general.
However, a phase transition of the C$_{\text{i}}$-dominated into a clearly C$_{\text{s}}$-dominated structure is observed.
The amount of substitutionally occupied C atoms increases with increasing temperature.
-Isolated structures of stretched SiC adjusted to the c-Si host with respect to the lattice constant and alignement are formed.
+Isolated structures of stretched SiC adjusted to the c-Si host with respect to the lattice constant and alignment are formed.
Entropic contributions are assumed to be responsible for these structures at elevated temperatures that deviate from the ground state at 0 K.
Results of the MD simulations at different temperatures and C concentrations can be correlated to experimental findings.
Some of the key properties are listed in Table~\ref{table:sic:properties} and compared to other technologically relevant semiconductor materials.
Despite the lower charge carrier mobilities for low electric fields SiC outperforms Si concerning all other properties.
The wide band gap, large breakdown field and high saturation drift velocity make SiC an ideal candidate for high-temperature, high-power and high-frequency electronic devices exhibiting high efficiency~\cite{wesch96,morkoc94,casady96,capano97,pensl93,park98,edgar92}.
-In addition the high thermal conductivity enables the implementation of small-sized electronic devices enduring increased power densites.
+In addition the high thermal conductivity enables the implementation of small-sized electronic devices enduring increased power densities.
Its formidable mechanical stability, heat resistant, radiation hardness and low neutron capture cross section allow operation in harsh and radiation-hard environments~\cite{capano97}.
Despite high-temperature operations the wide band gap also allows the use of SiC in optoelectronic devices.
The focus of SiC based applications, however, is in the area of solid state electronics experiencing revolutionary performance improvements enabled by its capabilities.
These devices include ultraviolet (UV) detectors, high power radio frequency (RF) amplifiers, rectifiers and switching transistors as well as microelectromechanical system (MEMS) applications.
-For UV dtectors the wide band gap is useful for realizing low photodiode dark currents as well as sensors that are blind to undesired near-infrared wavelenghts produced by heat and solar radiation.
+For UV detectors the wide band gap is useful for realizing low photodiode dark currents as well as sensors that are blind to undesired near-infrared wavelengths produced by heat and solar radiation.
These photodiodes serve as excellent sensors applicable in the monitoring and control of turbine engine combustion.
The low dark currents enable the use in X-ray, heavy ion and neutron detection in nuclear reactor monitoring and enhanced scientific studies of high-energy particle collisions as well as cosmic radiation.
The low neutron capture cross section and radiation hardness favors its use in detector applications.
For instance, SiC based solid state transmitters hold great promise for High Definition Television (HDTV) broadcast stations abandoning the reliance on tube-based technology for high-power transmitters significantly reducing the size of such transmitters and long-term maintenance costs.
The high breakdown field of SiC compared to Si allows the blocking voltage region of a device to be designed roughly 10 times thinner and 10 times heavier doped, resulting in a decrease of the blocking region resistance by a factor of 100 and a much faster switching behavior.
Thus, rectifier diodes and switching transistors with higher switching frequencies and much greater efficiencies can be realized and exploited in highly efficient power converters.
-Therefor, SiC constitutes a promising candidate to become the key technology towards an extensive development and use of regenerative energies and elctromobility.
+Therefor, SiC constitutes a promising candidate to become the key technology towards an extensive development and use of regenerative energies and electromobility.
Beside the mentioned electrical capabilities the mechanical stability, which is almost as hard as diamond, and chemical inertness almost suggest SiC to be used in MEMS designs.
Among the different polytypes of SiC, the cubic phase shows a high electron mobility and the highest break down field as well as saturation drift velocity.
\begin{center}
\includegraphics[width=0.35\columnwidth]{sic_unit_cell.eps}
\end{center}
-\caption[3C-SiC unit cell.]{3C-SiC unit cell. Yellow and grey spheres correpsond to Si and C atoms respectively. Covalent bonds are illustrated by blue lines.}
+\caption[3C-SiC unit cell.]{3C-SiC unit cell. Yellow and gray spheres correspond to Si and C atoms respectively. Covalent bonds are illustrated by blue lines.}
\label{fig:sic:unit_cell}
\end{figure}
Its unit cell is shown in Fig.~\ref{fig:sic:unit_cell}.
In the van Arkel apparatus \cite{arkel25}, Si and C containing gases like methylchlorosilanes \cite{moers31} and silicon tetrachloride \cite{kendall53} are pyrolitically decomposed and SiC is deposited on heated carbon rods in a vapor growth process.
Typical deposition temperatures are in the range between \unit[1400]{$^{\circ}$C} and \unit[1600]{$^{\circ}$C} while studies up to \unit[2500]{$^{\circ}$C} have been performed.
-The obtained polycrystalline material consists of small crystal grains with a size of several hunderd microns stated to be mainly of the cubic polytype.
+The obtained polycrystalline material consists of small crystal grains with a size of several hundred microns stated to be mainly of the cubic polytype.
A significant breakthrough was made in 1955 by Lely, who proposed a sublimation process for growing higher purity bulk SiC single crystals \cite{lely55}.
In the so called Lely process, a tube of porous graphite is surrounded by polycrystalline SiC as gained by previously described processes.
-Heating the hollow carbon cylinder to \unit[2500]{$^{\circ}$C} leads to sublimation of the material at the hot outer wall and diffusion through the porous graphite tube followed by an uncontrolled crystallization on the slightly cooler parts of the inner graphite cavity resulting in the formation of randomly sized, hexagonally shaped platelets, which exibit a layered structure of various alpha polytypes with equal \hkl{0001} orientation.
+Heating the hollow carbon cylinder to \unit[2500]{$^{\circ}$C} leads to sublimation of the material at the hot outer wall and diffusion through the porous graphite tube followed by an uncontrolled crystallization on the slightly cooler parts of the inner graphite cavity resulting in the formation of randomly sized, hexagonally shaped platelets, which exhibit a layered structure of various alpha polytypes with equal \hkl{0001} orientation.
Subsequent research \cite{tairov78,tairov81} resulted in the implementation of a seeded growth sublimation process wherein only one large crystal of a single polytype is grown.
In the so called modified Lely or modified sublimation process nucleation occurs on a SiC seed crystal located at the top or bottom of a cylindrical growth cavity.
The porous material constitutes a severe source of contamination, e.g. with the dopants N, B and Al, which is particularly effective at low temperatures due to the low growth rate.
Since the vapor pressure of Si is much higher than that of C, a careful manipulation of the Si vapor content above the seed crystal is required.
Additionally, to preserve epitaxial growth conditions, graphitization of the seed crystal has to be avoided.
-Avoiding defects constitutes a mojor difficulty.
+Avoiding defects constitutes a major difficulty.
These defects include growth spirals (stepped screw dislocations), subgrain boundaries and twins as well as micropipes (micron sized voids extending along the c axis of the crystal) and 3C inclusions at the seed crystal in hexagonal growth systems.
Micropipe-free growth of 6H-SiC has been realized by a reduction of the temperature gradient in the sublimation furnace resulting in near-equilibrium growth conditions in order to avoid stresses, which is, however, accompanied by a reduction of the growth rate \cite{schulze98}.
Further efforts have to be expended to find relations between the growth parameters, the kind of polytype and the occurrence and concentration of defects, which are of fundamental interest and might help to improve the purity of the bulk materials.
This results in the thermodynamically favored growth of a single phase due to the uni-directional contraction of Si-C-Si bond chains perpendicular to the terrace steps edges during carbonization and the fast growth parallel to the terrace edges during growth under Si rich conditions \cite{kitabatake97}.
By MBE, lower process temperatures than these typically employed in CVD have been realized \cite{hatayama95,henke95,fuyuki97,takaoka98}, which is essential for limiting thermal stresses and to avoid resulting substrate bending, a key issue in obtaining large area 3C-SiC surfaces.
In summary, the almost universal use of Si has allowed significant progress in the understanding of heteroepitaxial growth of SiC on Si.
-However, mismatches in the thermal expansion coefficient and the lattice parameter cause a considerably high concentration of various defects, which is responsible for structural and electrical qualities that are not yet statisfactory.
+However, mismatches in the thermal expansion coefficient and the lattice parameter cause a considerably high concentration of various defects, which is responsible for structural and electrical qualities that are not yet satisfactory.
The alternative attempt to grow SiC on SiC substrates has shown to drastically reduce the concentration of defects in deposited layers.
By CVD, both, the 3C \cite{kong88,powell90} as well as the 6H \cite{kong88_2,powell90_2} polytype could be successfully grown.
In order to obtain the homoepitaxially grown 6H polytype, off-axis 6H-SiC wafers are required as a substrate \cite{kimoto93}.
%In the so called step-controlled epitaxy, lateral growth proceeds from atomic steps without the necessity of preceding nucleation events.
-Investigations indicate that in the so-called step-controlled epitaxy, crystal growth proceeds through the adsorbtion of Si species at atomic steps and their carbonization by hydrocarbon molecules.
+Investigations indicate that in the so-called step-controlled epitaxy, crystal growth proceeds through the adsorption of Si species at atomic steps and their carbonization by hydrocarbon molecules.
This growth mechanism does not require two-dimensional nucleation.
Instead, crystal growth is governed by mass transport, i.e. the diffusion of reactants in a stagnant layer.
In contrast, layers of the 3C polytype are formed on exactly oriented \hkl(0 0 0 1) 6H-SiC substrates by two-dimensional nucleation on terraces.
These films show a high density of double positioning boundary (DPB) defects, which is a special type of twin boundary arising at the interface of regions that occupy one of the two possible orientations of the hexagonal stacking sequence, which are rotated by \unit[60]{$^{\circ}$} relative to each other, respectively.
However, lateral 3C-SiC growth was also observed on low tilt angle off-axis substrates originating from intentionally induced dislocations \cite{powell91}.
Additionally, 6H-SiC was observed on clean substrates even for a tilt angle as low as \unit[0.1]{$^{\circ}$} due to low surface mobilities that facilitate arriving molecules to reach surface steps.
-Thus, 3C nucleation is assumed as a result of migrating Si and C cointaining molecules interacting with surface disturbances by a yet unknown mechanism, in contrast to a model \cite{ueda90}, in which the competing 6H versus 3C growth depends on the density of surface steps.
+Thus, 3C nucleation is assumed as a result of migrating Si and C containing molecules interacting with surface disturbances by a yet unknown mechanism, in contrast to a model \cite{ueda90}, in which the competing 6H versus 3C growth depends on the density of surface steps.
Combining the fact of a well defined 3C lateral growth direction, i.e. the tilt direction, and an intentionally induced dislocation enables the controlled growth of a 3C-SiC film mostly free of DPBs \cite{powell91}.
Lower growth temperatures, a clean growth ambient, in situ control of the growth process, layer-by-layer deposition and the possibility to achieve dopant profiles within atomic dimensions due to the reduced diffusion at low growth temperatures reveal MBE as a promising technique to produce SiC epitaxial layers.
The beneficial epitaxial relation of substrate and film limits the structural difference between the two polytypes in two out of six layers with respect to the stacking sequence along the c axis.
Homoepitaxial growth of 3C-SiC by GSMBE was realized for the first time by atomic layer epitaxy (ALE) utilizing the periodical change in the surface superstructure by the alternating supply of the source gases, which determines the growth rate giving atomic level control in the growth process \cite{fuyuki89}.
The cleaned substrate surface shows a C terminated $(2\times 2)$ pattern at \unit[1000]{$^{\circ}$C}, which turns into a $(3\times 2)$ pattern when Si$_2$H$_6$ is introduced and it is maintained after the supply is stopped.
-A more detailed investigation showed the formation of a preceeding $(2\times 1)$ and $(5\times 2)$ pattern within the exposure to the Si containing gas \cite{yoshinobu90,fuyuki93}.
+A more detailed investigation showed the formation of a preceding $(2\times 1)$ and $(5\times 2)$ pattern within the exposure to the Si containing gas \cite{yoshinobu90,fuyuki93}.
The $(3\times 2)$ superstructure contains approximately 1.7 monolayers of Si atoms, crystallizing into 3C-SiC with a smooth and mirror-like surface after C$_2$H$_6$ is inserted accompanied by a reconstruction of the surface into the initial C terminated $(2\times 2)$ pattern.
A minimal growth rate of 2.3 monolayers per cycle exceeding the value of 1.7 is due to physically adsorbed Si atoms not contributing to the superstructure.
To realize single monolayer growth precise control of the gas supply to form the $(2\times 1)$ structure is required.
However, accurate layer-by-layer growth is achieved under certain conditions, which facilitate the spontaneous desorption of an additional layer of one atom species by supply of the other species \cite{hara93}.
Homoepitaxial growth of the 6H polytype has been realized on off-oriented substrates utilizing simultaneous supply of the source gases \cite{tanaka94}.
Depending on the gas flow ratio either 3C island formation or step flow growth of the 6H polytype occurs, which is explained by a model including aspects of enhanced surface mobilities of adatoms on a $(3\times 3)$ reconstructed surface.
-Due to the strong adsorption of atomic hydrogen \cite{allendorf91} decomposited of the gas phase reactants at low temperatures, however, there seems to be no benefit of GSMBE compared to CVD.
+Due to the strong adsorption of atomic hydrogen \cite{allendorf91} decomposed of the gas phase reactants at low temperatures, however, there seems to be no benefit of GSMBE compared to CVD.
Next to lattice imperfections, incorporated hydrogen effects the surface mobility of the adsorbed species \cite{eaglesham93} setting a minimum limit for the growth temperature, which would preferably be further decreased in order to obtain sharp doping profiles.
Thus, growth rates must be adjusted to be lower than the desorption rate of hydrogen, which leads to very low deposition rates at low temperatures.
-SSMBE, by supplying the atomic species to be deposited by evaporation of a solid, presumably constitutes the preffered method in order to avoid the problems mentioned above.
+SSMBE, by supplying the atomic species to be deposited by evaporation of a solid, presumably constitutes the preferred method in order to avoid the problems mentioned above.
Although, in the first experiments, temperatures still above \unit[1100]{$^{\circ}$C} were necessary to epitaxially grow 3C-SiC films on 6H-SiC substrates \cite{kaneda87}, subsequent attempts succeeded in growing mixtures of twinned 3C-SiC and 6H-SiC films on off-axis \hkl(0001) 6H-SiC wafers at temperatures between \unit[800]{$^{\circ}$C} and \unit[1000]{$^{\circ}$C} \cite{fissel95,fissel95_apl}.
In the latter approach, as in GSMBE, excess Si atoms, which are controlled by the Si/C flux ratio, result in the formation of a Si adlayer and the formation of a non-stoichiometric, reconstructed surface superstructure, which influences the mobility of adatoms and, thus, has a decisive influence on the growth mode, polytype and crystallinity \cite{fissel95,fissel96,righi03}.
Therefore, carefully controlling the Si/C ratio could be exploited to obtain definite heterostructures of different SiC polytypes providing the possibility for band gap engineering in SiC materials.
To summarize, much progress has been achieved in SiC thin film growth during the last few years.
-However, the frequent occurence of defects such as dislocations, twins and double positioning boundaries limit the structural and electrical quality of large SiC films.
+However, the frequent occurrence of defects such as dislocations, twins and double positioning boundaries limit the structural and electrical quality of large SiC films.
Solving this issue remains a challenging problem necessary to drive SiC for potential applications in high-performance electronic device production \cite{wesch96}.
\subsection{Ion beam synthesis of cubic silicon carbide}
\label{subsection:ibs}
-Although tremendous progress has been achieved in the above-mentioned growth methods during the last decades, available wafer dimensions and crystal qualities are not yet statisfactory.
+Although tremendous progress has been achieved in the above-mentioned growth methods during the last decades, available wafer dimensions and crystal qualities are not yet satisfactory.
Thus, alternative approaches to fabricate SiC have been explored.
The ion beam synthesis (IBS) technique, i.e. high-dose ion implantation followed by a high-temperature annealing step, turned out to constitute a promising method to directly form compound layers of high purity and accurately controllable depth and stoichiometry.
A short chronological summary of the IBS of SiC and its origins is presented in the following.
High-dose carbon implantation into crystalline silicon (c-Si) with subsequent or in situ annealing was found to result in SiC microcrystallites in Si \cite{borders71}.
-Rutherford backscattering spectrometry (RBS) and infrared (IR) spectroscopy investigations indicate a \unit[10]{at.\%} C concentration peak and the occurence of disordered C-Si bonds after implantation at room temperature (RT) followed by crystallization into SiC precipitates upon annealing demonstrated by a shift in the IR absorption band and the disappearance of the C profile peak in RBS.
+Rutherford backscattering spectrometry (RBS) and infrared (IR) spectroscopy investigations indicate a \unit[10]{at.\%} C concentration peak and the occurrence of disordered C-Si bonds after implantation at room temperature (RT) followed by crystallization into SiC precipitates upon annealing demonstrated by a shift in the IR absorption band and the disappearance of the C profile peak in RBS.
Implantations at different temperatures revealed a strong influence of the implantation temperature on the compound structure \cite{edelman76}.
Temperatures below \unit[500]{$^{\circ}$C} result in amorphous layers, which is transformed into polycrystalline 3C-SiC after \unit[850]{$^{\circ}$C} annealing.
Otherwise single crystalline 3C-SiC is observed for temperatures above \unit[600]{$^{\circ}$C}.
Overstoichiometric doses result in the formation of clusters of C, which do not contribute to SiC formation during annealing up to \unit[1200]{$^{\circ}$C} \cite{kimura82}.
The amount of formed SiC, however, increases with increasing implantation temperature.
The authors, thus, concluded that implantations at elevated temperatures lead to a reduction in the annealing temperatures required for the synthesis of homogeneous layers of SiC.
-In a comparative study of O, N and C implantation into Si, the absence of the formation of a stoichiometric SiC compound layer involving the transition of a Gaussian into a box-like C profile with respect to the implantation depth for the superstoichiometric C implantation and an annealing temeprature of \unit[1200]{$^{\circ}$C} in contrast to the O and N implantations, which successfully form homogeneous layers, has been observed \cite{reeson86}.
-This was attrubuted to the difference in the enthalpy of formation of the respective compound and the different mobility of the respective impurity in bulk Si.
+In a comparative study of O, N and C implantation into Si, the absence of the formation of a stoichiometric SiC compound layer involving the transition of a Gaussian into a box-like C profile with respect to the implantation depth for the superstoichiometric C implantation and an annealing temperature of \unit[1200]{$^{\circ}$C} in contrast to the O and N implantations, which successfully form homogeneous layers, has been observed \cite{reeson86}.
+This was attributed to the difference in the enthalpy of formation of the respective compound and the different mobility of the respective impurity in bulk Si.
Thus, higher annealing temperatures and longer annealing times were considered necessary for the formation of homogeneous SiC layers.
Indeed, for the first time, buried homogeneous and stoichiometric epitaxial 3C-SiC layers embedded in single crystalline Si were obtained by the same group consequently applying annealing temperatures of \unit[1405]{$^{\circ}$C} for \unit[90]{min} and implantation temperatures of approximately \unit[550]{$^{\circ}$C} \cite{reeson87}.
The necessity of the applied extreme temperature and time scale is attributed to the stability of substitutional C within the Si matrix being responsible for high activation energies necessary to dissolve such precipitates and, thus, allow for redistribution of the implanted C atoms.
Hence, to obtain sharp interfaces and single-crystalline SiC layers temperatures between \unit[400]{$^{\circ}$C} and \unit[600]{$^{\circ}$C} have to be used.
Indeed, reasonable results were obtained at \unit[500]{$^{\circ}$C} \cite{lindner98} and even better interfaces were observed for \unit[450]{$^{\circ}$C} \cite{lindner99_2}.
To further improve the interface quality and crystallinity a two-temperature implantation technique was developed \cite{lindner99}.
-To form a narrow, box-like density profile of oriented SiC nanocrystals \unit[93]{\%} of the total dose of \unit[$8.5\cdot 10^{17}$]{cm$^{-2}$} is implanted at \unit[500]{$^{\circ}$C}.
-The remaining dose is implanted at \unit[250]{$^{\circ}$C}, which leads to the formation of amorphous zones above and below the SiC precipitate layer and the desctruction of SiC nanocrystals within these zones.
+To form a narrow, box-like density profile of oriented SiC nanocrystals, \unit[93]{\%} of the total dose of \unit[$8.5\cdot 10^{17}$]{cm$^{-2}$} is implanted at \unit[500]{$^{\circ}$C}.
+The remaining dose is implanted at \unit[250]{$^{\circ}$C}, which leads to the formation of amorphous zones above and below the SiC precipitate layer and the destruction of SiC nanocrystals within these zones.
After annealing for \unit[10]{h} at \unit[1250]{$^{\circ}$C} a homogeneous, stoichiometric SiC layer with sharp interfaces is formed.
Fig. \ref{fig:sic:hrem_sharp} shows the respective high resolution transmission electron microscopy micrographs.
\begin{figure}[t]
Similar to phosphorous and boron, which exclusively use self-interstitials as a diffusion vehicle, the diffusion of C atoms is expected to obey the same mechanism.
Indeed, enhanced C diffusion was observed in the presence of self-interstitial supersaturation \cite{kalejs84} indicating an appreciable diffusion component involving self-interstitials and only a negligible contribution by vacancies.
Substitutional C and interstitial Si react into a C-Si complex forming a dumbbell structure oriented along a crystallographic \hkl<1 0 0> direction on a regular Si lattice site.
-This structure, the so called C-Si \hkl<1 0 0> dumbbell structure, was initially suspected by local vibrational mode absorption \cite{bean70} and finally verified by electron paramegnetic resonance \cite{watkins76} studies on irradiated Si substrates at low temperatures.
+This structure, the so called C-Si \hkl<1 0 0> dumbbell structure, was initially suspected by local vibrational mode absorption \cite{bean70} and finally verified by electron paramagnetic resonance \cite{watkins76} studies on irradiated Si substrates at low temperatures.
Measuring the annealing rate of the defect as a function of temperature reveals barriers for migration ranging from \unit[0.73]{eV} \cite{song90} to \unit[0.87]{eV} \cite{tipping87}, which is highly mobile compared to substitutional C.
% diffusion pathway?
% expansion of the lattice (positive strain)
% lattice location of implanted carbon
Radiation damage introduced during implantation and a high concentration of the implanted species, which results in the reduction of the topological constraint of the host lattice imposed on the implanted species, can affect the manner of impurity incorporation.
The probability of finding C, which will be most stable at sites for which the number of neighbors equals the natural valence, i.e. substitutionally on a regular Si site of a perfect lattice, is, thus, reduced at substitutional lattice sites and likewise increased at interstitial sites.
-Indeed, x-ray rocking curves reveal a positive lattice strain, which is decreased but still remains with increasing annealing temeprature, indicating the location of the majority of implanted C atoms at interstitial sites \cite{isomae93}.
+Indeed, x-ray rocking curves reveal a positive lattice strain, which is decreased but still remains with increasing annealing temperature, indicating the location of the majority of implanted C atoms at interstitial sites \cite{isomae93}.
Due to the absence of dislocations in the implanted region interstitial C is assumed to prevent clustering of implantation-induced Si self-interstitials by agglomeration of C-Si interstitials or the formation of SiC precipitates accompanied by a relaxation of the lattice strain.
% link to strain engineering
However, there is great interest to incorporate C onto substitutional lattice sites, which results in a contraction of the Si lattice due to the smaller covalent radius of C compared to Si \cite{baker68}, causing tensile strain, which is applied to the Si lattice.
Thus, substitutional C enables strain engineering of Si and Si/Si$_{1-x}$Ge$_x$ heterostructures \cite{yagi02,chang05,kissinger94,osten97}, which is used to increase charge carrier mobilities in Si as well as to adjust its band structure \cite{soref91,kasper91}.
% increase of C at substitutional sites
-Epitaxial layers with \unit[1.4]{at.\%} of substitutional C have been successfully synthesized in preamorphized Si$_{0.86}$Ge$_{0.14}$ layers, which were grown by CVD on Si substrates, using multiple-energy C implantation followed by solid-physe epitaxial regrowth at \unit[700]{$^{\circ}$C} \cite{strane93}.
+Epitaxial layers with \unit[1.4]{at.\%} of substitutional C have been successfully synthesized in preamorphized Si$_{0.86}$Ge$_{0.14}$ layers, which were grown by CVD on Si substrates, using multiple-energy C implantation followed by solid-phase epitaxial regrowth at \unit[700]{$^{\circ}$C} \cite{strane93}.
The tensile strain induced by the C atoms is found to compensate the compressive strain present due to the Ge atoms.
Studies on the thermal stability of Si$_{1-y}$C$_y$/Si heterostructures formed in the same way and equal C concentrations showed a loss of substitutional C accompanied by strain relaxation for temperatures ranging from \unit[810-925]{$^{\circ}$C} and the formation of spherical 3C-SiC precipitates with diameters of \unit[2-4]{nm}, which are incoherent but aligned to the Si host \cite{strane94}.
During the initial stages of precipitation C-rich clusters are assumed, which maintain coherency with the Si matrix and the associated biaxial strain.
-Using this technique a metastable solubility limit was achieved, which corresponds to a C concentration exceeding the solid solubility limit at the Si melting point by nearly three orders of magnitude and, furthermore, a reduction of the defect denisty near the metastable solubility limit is assumed if the regrowth temperature is increased by rapid thermal annealing \cite{strane96}.
+Using this technique a metastable solubility limit was achieved, which corresponds to a C concentration exceeding the solid solubility limit at the Si melting point by nearly three orders of magnitude and, furthermore, a reduction of the defect density near the metastable solubility limit is assumed if the regrowth temperature is increased by rapid thermal annealing \cite{strane96}.
Since high temperatures used in the solid-phase epitaxial regrowth method promotes SiC precipitation, other groups realized substitutional C incorporation for strained Si$_{1-y}$C$_y$/Si heterostructures \cite{iyer92,fischer95,powell93,osten96,osten99,laveant2002} or partially to fully strain-compensated (even inversely distorted \cite{osten94_2}) Si$_{1-x-y}$Ge$_x$C${_y}$ layers on Si \cite{eberl92,powell93_2,osten94,dietrich94} by MBE.
Investigations reveal a strong dependence of the growth temperature on the amount of substitutionally incorporated C, which is increased for decreasing temperature accompanied by deterioration of the crystal quality \cite{osten96,osten99}.
While not being compatible to very-large-scale integration technology, C concentrations of \unit[2]{\%} and more have been realized \cite{laveant2002}.
\label{section:assumed_prec}
Although high-quality films of single-crystalline 3C-SiC can be produced by means of IBS the precipitation mechanism in bulk Si is not yet fully understood.
-Indeed, closely investigating the large amount of literature pulled up in the last two sections and a cautios combination of some of the findings reveals controversial ideas of SiC formation, which are reviewed in more detail in the following.
+Indeed, closely investigating the large amount of literature pulled up in the last two sections and a cautious combination of some of the findings reveals controversial ideas of SiC formation, which are reviewed in more detail in the following.
High resolution transmission electron microscopy (HREM) investigations of C-implanted Si at room temperature followed by rapid thermal annealing (RTA) show the formation of C-Si dumbbell agglomerates, which are stable up to annealing temperatures of about \unit[700-800]{$^{\circ}$C}, and a transformation into 3C-SiC precipitates at higher temperatures \cite{werner96,werner97}.
-The precipitates with diamateres between \unit[2]{nm} and \unit[5]{nm} are incorporated in the Si matrix without any remarkable strain fields, which is explained by the nearly equal atomic density of C-Si agglomerates and the SiC unit cell.
+The precipitates with diameters between \unit[2]{nm} and \unit[5]{nm} are incorporated in the Si matrix without any remarkable strain fields, which is explained by the nearly equal atomic density of C-Si agglomerates and the SiC unit cell.
Implantations at \unit[500]{$^{\circ}$C} likewise suggest an initial formation of C-Si dumbbells on regular Si lattice sites, which agglomerate into large clusters \cite{lindner99_2}.
The agglomerates of such dimers, which do not generate lattice strain but lead to a local increase of the lattice potential \cite{werner96,werner97}, are indicated by dark contrasts and otherwise undisturbed Si lattice fringes in HREM, as can be seen in Fig.~\ref{fig:sic:hrem:c-si}.
\begin{figure}[t]
The insignificantly lower Si density of SiC of approximately \unit[3]{\%} compared to c-Si results in the emission of only a few excess Si atoms.
The same mechanism was identified by high resolution x-ray diffraction \cite{eichhorn99}.
For implantation temperatures of \unit[500]{$^{\circ}$C} C-Si dumbbells agglomerate in an initial stage followed by the additional appearance of aligned SiC precipitates in a slightly expanded Si region with increasing dose.
-The precipitation mechanism based on a preceeding dumbbell agglomeration as indicated by the above-mentioned experiemnts is schematically displayed in Fig.~\ref{fig:sic:db_agglom}.
+The precipitation mechanism based on a preceding dumbbell agglomeration as indicated by the above-mentioned experiments is schematically displayed in Fig.~\ref{fig:sic:db_agglom}.
\begin{figure}[t]
\begin{center}
\subfigure[]{\label{fig:sic:db_agglom:seq01}\includegraphics[width=0.30\columnwidth]{sic_prec_seq_01.eps}}
\chapter{Details of simulation parameters and test calculations}
\label{chapter:simulation}
-All calculations are carried out utilizing the supercell approach, which means that the simulation cell contains a multiple of unti cells and periodic boundary conditions are imposed on the boundaries of that simulation cell.
+All calculations are carried out utilizing the supercell approach, which means that the simulation cell contains a multiple of unit cells and periodic boundary conditions are imposed on the boundaries of that simulation cell.
Strictly, these supercells become the unit cells, which, by a periodic sequence, compose the bulk material that is actually investigated by this approach.
Thus, importance need to be attached to the construction of the supercell.
Three basic types of supercells to compose the initial Si bulk lattice, which can be scaled by integers in the different directions, are considered.
An energy cut-off of \unit[300]{eV} is used to expand the wave functions into the plane-wave basis.
Sampling of the Brillouin zone is restricted to the $\Gamma$ point.
Spin polarization has been fully accounted for.
-The electronic ground state is calculated by an interative Davidson scheme \cite{davidson75} until the difference in total energy of two subsequent iterations is below \unit[$10^{-4}$]{eV}.
+The electronic ground state is calculated by an iterative Davidson scheme \cite{davidson75} until the difference in total energy of two subsequent iterations is below \unit[$10^{-4}$]{eV}.
Defect structures and the migration paths have been modeled in cubic supercells of type 3 containing 216 Si atoms.
-The conjugate gradiant algorithm is used for ionic relaxation.
+The conjugate gradient algorithm is used for ionic relaxation.
Migration paths are determined by the modified version of the CRT method as explained in section \ref{section:basics:migration}.
-The cell volume and shape is allowed to change using the pressure control algorithm of Parinello and Rahman \cite{parrinello81} in order to realize constant pressure simulations.
+The cell volume and shape is allowed to change using the pressure control algorithm of Parrinello and Rahman \cite{parrinello81} in order to realize constant pressure simulations.
Due to restrictions by the {\textsc vasp} code, {\em ab initio} MD could only be performed at constant volume.
-In MD simulations the equations of motion are integrated by a fourth order predictor corrector algorithm for a timestep of \unit[1]{fs}.
+In MD simulations the equations of motion are integrated by a fourth order predictor corrector algorithm for a time step of \unit[1]{fs}.
% todo - point defects are calculated for the neutral charge state.
\hline
\end{tabular}
\end{center}
-\caption[Equilibrium lattice constants and cohesive energies of fully relaxed structures of Si, C (diamond) and 3C-SiC for different potentials and XC functionals.]{Equilibrium lattice constants and cohesive energies of fully relaxed structures of Si, C (diamond) and 3C-SiC for different potentials (ultr-soft PP and PAW) and XC functionals (LDA and GGA). Deviations of the respective values from experimental values are given. Values are in good (green), fair (orange) and poor (red) agreement.}
+\caption[Equilibrium lattice constants and cohesive energies of fully relaxed structures of Si, C (diamond) and 3C-SiC for different potentials and XC functionals.]{Equilibrium lattice constants and cohesive energies of fully relaxed structures of Si, C (diamond) and 3C-SiC for different potentials (ultra-soft PP and PAW) and XC functionals (LDA and GGA). Deviations of the respective values from experimental values are given. Values are in good (green), fair (orange) and poor (red) agreement.}
\label{table:simulation:potxc}
\end{table}
Table \ref{table:simulation:potxc} shows the lattice constants and cohesive energies obtained for the fully relaxed structures with respect to the utilized potential and XC functional.
The 3C-SiC calculations employing the PAW method in conjunction with the LDA suffers from the general problem inherent to LDA, i.e. overestimated binding energies.
Thus, the PAW \& LDA combination is not pursued.
Since the lattice constant and cohesive energy of 3C-SiC calculated by the PAW method using the GGA are not improved compared to the ultra-soft pseudopotential calculations using the same XC functional, this concept is likewise stopped.
-To conclude, the combination of ultr-soft pseudopotentials and the GGA XC functional are considered the optimal choice for the present study.
+To conclude, the combination of ultra-soft pseudopotentials and the GGA XC functional are considered the optimal choice for the present study.
\subsection{Lattice constants and cohesive energies}
\end{table}
Table \ref{table:simulation:paramf} shows the respective results and deviations from experiment.
A nice agreement with experimental results is achieved.
-Clearly, a competent parameter set is found, which is capabale of describing the C/Si system by {\em ab initio} calculations.
+Clearly, a competent parameter set is found, which is capable of describing the C/Si system by {\em ab initio} calculations.
% todo - ref for experimental values!
\label{section:classpotmd}
The classical potential MD method is much less computationally costly compared to the highly accurate quantum-mechanical method.
-Thus, the method is capable of performing structural optimizations on large systems and MD calulations may be used to model a system over long time scales.
+Thus, the method is capable of performing structural optimizations on large systems and MD calculations may be used to model a system over long time scales.
Defect structures are modeled in a cubic supercell (type 3) of nine Si lattice constants in each direction containing 5832 Si atoms.
Reproducing the SiC precipitation was attempted in cubic c-Si supercells, which have a size up to 31 Si unit cells in each direction consisting of 238328 Si atoms.
A Tersoff-like bond order potential by Erhart and Albe (EA) \cite{albe_sic_pot} is used to describe the atomic interaction.
Integration of the equations of motion is realized by the velocity Verlet algorithm \cite{verlet67} using a fixed time step of \unit[1]{fs}.
For structural relaxation of defect structures the same algorithm is utilized with the temperature set to zero Kelvin.
This also applies for the relaxation of structures within the CRT calculations to find migration pathways.
-In the latter case the time constant of the Berendsen thermostat is set to \unit[1]{fs} in order to achieve direct velocity scaling, which corresponds to a steepest descent minimazation driving the system into a local minimum, if the temperature is set to zero Kelvin.
+In the latter case the time constant of the Berendsen thermostat is set to \unit[1]{fs} in order to achieve direct velocity scaling, which corresponds to a steepest descent minimization driving the system into a local minimum, if the temperature is set to zero Kelvin.
However, in some cases a time constant of \unit[100]{fs} turned out to result in lower barriers.
Defect structures as well as the simulations modeling the SiC precipitation are performed in the isothermal-isobaric $NpT$ ensemble.
\label{fig:simulation:pc_sic-prec}
\end{figure}
Fig. \ref{fig:simulation:pc_sic-prec} shows the radial distribution of the obtained precipitate configuration.
-The Si-Si radial distribution for both, plain c-Si and the precipitate configuration show a maximum at a distance of \unit[0.235]{nm}, which is the distance of next neighboured Si atoms in c-Si.
+The Si-Si radial distribution for both, plain c-Si and the precipitate configuration show a maximum at a distance of \unit[0.235]{nm}, which is the distance of next neighbored Si atoms in c-Si.
Although no significant change of the lattice constant of the surrounding c-Si matrix was assumed, surprisingly, there is no change at all within observational accuracy.
Looking closer at higher order Si-Si peaks might even allow the guess of a slight increase of the lattice constant compared to the plain c-Si structure.
A new Si-Si peak arises at \unit[0.307]{nm}, which is identical to the peak of the C-C distribution around that value.
The bumps of the Si-Si distribution at higher distances marked by the green arrows can be explained in the same manner.
They correspond to the fourth and sixth next neighbor distance in 3C-SiC.
It is easily identifiable how these C-C peaks, which imply Si pairs at same distances inside the precipitate, contribute to the bumps observed in the Si-Si distribution.
-The Si-Si and C-C peak at \unit[0.307]{nm} enables the determination of the lattic constant of the embedded 3C-SiC precipitate.
+The Si-Si and C-C peak at \unit[0.307]{nm} enables the determination of the lattice constant of the embedded 3C-SiC precipitate.
A lattice constant of \unit[4.34]{\AA} compared to \unit[4.36]{\AA} for bulk 3C-SiC is obtained.
This is in accordance with the peak of Si-C pairs at a distance of \unit[0.188]{nm}.
Thus, the precipitate structure is slightly compressed compared to the bulk phase.
\end{equation}
with the notation used in Table \ref{table:simulation:sic_prec}.
Here, $a_{\text{c-Si prec}}$ denotes the lattice constant of the surrounding crystalline Si and $a_{\text{3C-SiC prec}}$ is the lattice constant of the precipitate.
-The lattice constant of plain c-Si at \unit[20]{$^{\circ}$C} can be determined more accurately by the side lengthes of the simulation box of an equilibrated structure instead of using the radial distribution data.
+The lattice constant of plain c-Si at \unit[20]{$^{\circ}$C} can be determined more accurately by the side lengths of the simulation box of an equilibrated structure instead of using the radial distribution data.
By this, a value of $a_{\text{plain c-Si}}=5.439\,\text{\AA}$ is obtained.
The same lattice constant is assumed for the c-Si surrounding in the precipitate configuration $a_{\text{c-Si prec}}$ since peaks in the radial distribution match the ones of plain c-Si.
Using $a_{\text{3C-SiC prec}}=4.34\,\text{\AA}$ as observed from the radial distribution finally results in an increase of the initial volume by \unit[0.12]{\%}.
where $E$ is the total energy of the precipitate configuration at zero temperature.
An interfacial energy of \unit[2267.28]{eV} is obtained.
The amount of C atoms together with the observed lattice constant of the precipitate leads to a precipitate radius of \unit[29.93]{\AA}.
-Thus, the interface tension, given by the energy of the interface devided by the surface area of the precipitate is \unit[20.15]{eV/nm$^2$} or \unit[$3.23\times 10^{-4}$]{J/cm$^2$}.
-This value perfectly fits within the eperimentally estimated range of \unit[$2-8\times10^{-4}$]{J/cm$^2$} \cite{taylor93}.
+Thus, the interface tension, given by the energy of the interface divided by the surface area of the precipitate is \unit[20.15]{eV/nm$^2$} or \unit[$3.23\times 10^{-4}$]{J/cm$^2$}.
+This value perfectly fits within the experimentally estimated range of \unit[$2-8\times10^{-4}$]{J/cm$^2$} \cite{taylor93}.
Thus, the EA potential is considered an appropriate choice for the current study concerning the accurate description of the energetics of interfaces.
Furthermore, since the calculated interfacial energy is located in the lower part of the experimental range, the obtained interface structure might resemble an authentic configuration of an energetically favorable interface structure of a 3C-SiC precipitate in c-Si.
\subsubsection{Stability of the precipitate}
-To investigate the stability of the precipiate, which is assumed to be stable even at temperatures above the Si melting temperature, the configuration is heated up beyond the critical point, at which the Si melting transition occurs.
+To investigate the stability of the precipitate, which is assumed to be stable even at temperatures above the Si melting temperature, the configuration is heated up beyond the critical point, at which the Si melting transition occurs.
For this, the transition point of c-Si needs to be evaluated first.
According to the authors of the potential, the Si melting point is \degk{2450}.
-However, melting is not predicted to occur instantly after exceeding the melting point due to the additionally required transition enthalpy and hysteresis behaviour.
+However, melting is not predicted to occur instantly after exceeding the melting point due to the additionally required transition enthalpy and hysteresis behavior.
To determine the transition point, c-Si is heated up using a heating rate of \unit[1]{$^{\circ}$C/ps}.
\begin{figure}[tp]
\begin{center}
\label{fig:simulation:fe_and_t}
\end{figure}
Fig.~\ref{fig:simulation:fe_and_t} shows the total energy and temperature evolution in the region around the transition temperature.
-Indeed, a transition and the accompanied critical behaviour of the total energy is first observed at approximately \degk{3125}, which corresponds to \unit[128]{\%} of the Si melting temperature.
+Indeed, a transition and the accompanied critical behavior of the total energy is first observed at approximately \degk{3125}, which corresponds to \unit[128]{\%} of the Si melting temperature.
The difference in total energy is \unit[0.58]{eV} per atom corresponding to \unit[55.7]{kJ/mole}, which compares quite well to the Si enthalpy of melting of \unit[50.2]{kJ/mole}.
The precipitate structure is heated up using the same heating rate.
}
These propose the precipitation of SiC by agglomeration of \ci{} DBs followed by a sudden formation of SiC and otherwise a formation by successive accumulation of \cs{} via intermediate stretched SiC structures, which are coherent to the Si lattice.
To solve this controversy and contribute to the understanding of SiC precipitation in c-Si, a series of atomistic simulations is carried out.
-In the first part, intrinsic and C related point defects in c-Si as well as some selected diffusion processes of the C defect are investigated by means of first-principles quatum-mechanical calculations based on DFT and classical potential calculations employing a Tersoff-like analytical bond order potential.
+In the first part, intrinsic and C related point defects in c-Si as well as some selected diffusion processes of the C defect are investigated by means of first-principles quantum-mechanical calculations based on DFT and classical potential calculations employing a Tersoff-like analytical bond order potential.
Shortcomings of the computationally efficient though less accurate classical potential approach compared to the quantum-mechanical treatment are revealed.
The study proceeds investigating combinations of defect structures and related diffusion processes exclusively by the first-principles method.
The applicability of the utilized bond order potential for subsequent MD simulations is discussed.
Defects of C in Si are well described by both methods.
The \ci{} \hkl<1 0 0> DB is found to constitute the most favorable interstitial configuration in agreement with several theoretical \cite{burnard93,leary97,dal_pino93,capaz94,jones04} and experimental \cite{watkins76,song90} investigations.
Almost equal formation energies are predicted by the EA and DFT calculations for this defect.
-A small discrepancy in the resulting equilibrium structure with respect to the DFT method exists due to missing quantum-mechanical effects within the calssical treatment.
+A small discrepancy in the resulting equilibrium structure with respect to the DFT method exists due to missing quantum-mechanical effects within the classical treatment.
The high formation energies of the tetrahedral and hexagonal configuration obtained by classical potentials act in concert with the fact that these configurations are found unstable by the first-principles description.
The BC configuration turns out to be unstable relaxing into the \ci{} \hkl<1 1 0> DB configuration within the EA approach.
-This is supported by another {\em ab inito} study \cite{capaz94}, which in turn finds the BC configuration to be an intermediate saddle point structure of a possible migration path, which is \unit[2.1]{eV} higher than the \ci{} \hkl<1 0 0> DB structure.
+This is supported by another {\em ab initio} study \cite{capaz94}, which in turn finds the BC configuration to be an intermediate saddle point structure of a possible migration path, which is \unit[2.1]{eV} higher than the \ci{} \hkl<1 0 0> DB structure.
By quantum-mechanical calculations performed in this work, however, it turns out that the BC configuration constitutes a real local minimum if the electron spin is fully accounted for.
Indeed, spin polarization is absolutely necessary for the BC configuration resulting in a net magnetization of two electrons accompanied by a reduction of the total energy by \unit[0.3]{eV}.
The resulting spin up density is localized in a torus around the C perpendicular to the linear Si-C-Si bond.
Based on the above findings, it is concluded that modeling of the SiC precipitation by the EA potential might lead to trustable results.
Quantum-mechanical investigations of the mobility of the \ci{} \hkl<1 0 0> DB yield a migration barrier of \unit[0.9]{eV}, which excellently agrees to experimental values ranging from \unit[0.70]{eV} to \unit[0.87]{eV} \cite{lindner06,song90,tipping87}.
-The respective path correpsonds to a \ci{} \hkl[0 0 -1] DB migrating towards the next neighbored Si atom located in \hkl[1 1 -1] direction forming a \ci{} \hkl[0 -1 0] DB.
+The respective path corresponds to a \ci{} \hkl[0 0 -1] DB migrating towards the next neighbored Si atom located in \hkl[1 1 -1] direction forming a \ci{} \hkl[0 -1 0] DB.
The identified migration path involves a change in orientation of the DB.
-Thus, the same path explains the experimentally determined activation energies for reorientation of the DB ranging from \unit[0.77]{eV} \cite{watkins76} upto \unit[0.88]{eV} \cite{song90}.
+Thus, the same path explains the experimentally determined activation energies for reorientation of the DB ranging from \unit[0.77]{eV} \cite{watkins76} up to \unit[0.88]{eV} \cite{song90}.
Investigations based on the EA bond order potential suggest a migration involving an intermediate \ci{} \hkl<1 1 0> DB configuration.
Although different, starting and final configuration as well as the change in orientation of the \hkl<1 0 0> DB are equal to the identified pathway by the {\em ab initio} calculations.
However, barrier heights, which are overestimated by a factor of 2.4 to 3.5 depending on the character of migration, i.e. a single step or two step process, compared to the DFT results, are obtained.
Results of combinations of \ci{} and \cs{} revealed two additional metastable structures different to these obtained by a naive relaxation.
Small displacements result in a structure of a \hkl<1 1 0> C-C DB and in a structure of a twofold coordinated Si atom located in between two substitutional C atoms residing on regular Si lattice sites.
-Both structures are lower in energy compared to the respetive counterparts.
+Both structures are lower in energy compared to the respective counterparts.
These results, for the most part, compare well with results gained in previous studies \cite{leary97,capaz98,liu02} and show an astonishingly good agreement with experiment \cite{song90_2}.
Again, spin polarized calculations are revealed necessary.
A net magnetization of two electrons is observed for the \hkl<1 1 0> C-C DB configuration, which constitutes the ground state.
Limitations of the MD technique in addition to overestimated bond strengths due to the short range potential are identified to be responsible.
The approach of using increased temperatures during C insertion is followed to work around this problem termed {\em potential enhanced slow phase space propagation}.
-Higher temperatures are justified for severeal reasons.
-Elevated temperatures are expected to compensate the overestimated diffusion barriers and accelerate strcutural evolution.
-In addition, formation of SiC is also observed at higher implantation temperatures \cite{nejim95,lindner01} and temperatures in the implantation region is definetly higher than the temperature determined experimentally at the surface of the sample.
+Higher temperatures are justified for several reasons.
+Elevated temperatures are expected to compensate the overestimated diffusion barriers and accelerate structural evolution.
+In addition, formation of SiC is also observed at higher implantation temperatures \cite{nejim95,lindner01} and temperatures in the implantation region is definitely higher than the temperature determined experimentally at the surface of the sample.
Furthermore, the present study focuses on structural transitions in a system far from equilibrium.
No significant change is observed for high C concentrations at increased temperatures.
The amorphous phase is maintained.
The huge amount of damage hampers identification of investigated structures, which in many cases lost the alignment to the c-Si host.
-Obviously, inccorporation of a high quantity of C into a small volume within a short period of time creates damage, which decelerates structural evolution.
+Obviously, incorporation of a high quantity of C into a small volume within a short period of time creates damage, which decelerates structural evolution.
For the low C concentrations, time scales are still too low to observe C agglomeration.
However, a phase transition of the C$_{\text{i}}$-dominated into a clearly C$_{\text{s}}$-dominated structure is observed.
The amount of \cs{} increases with increasing temperature.
Initially, quantum-mechanical investigations suggest agglomeration of \ci{} defects that form energetically favorable configurations by an effective stress compensation.
Low barriers of migration are found except for transitions into the ground-state configuration, which is composed of a strong C-C bond.
-Thus, agglomeration of \ci{} in the absence of C clsutering is expected.
-These initial results suggest a conversion mechansim based on the agglomeration of \ci{} defects followed by a sudden precipitation once a critical size is reached.
-However, subsequent investigations of structures that are particularly conceivable under conditions prevalant in IBS and at elevated temperatures show \cs{} to occur in all probability.
+Thus, agglomeration of \ci{} in the absence of C clustering is expected.
+These initial results suggest a conversion mechanism based on the agglomeration of \ci{} defects followed by a sudden precipitation once a critical size is reached.
+However, subsequent investigations of structures that are particularly conceivable under conditions prevalent in IBS and at elevated temperatures show \cs{} to occur in all probability.
The transition from the ground state of a single C atom incorporated into otherwise perfect c-Si, i.e. the \ci{} \hkl<1 0 0> DB, into a configuration of \cs{} next to a \si{} atom exhibits an activation energy lower than the one for the diffusion of the highly mobile \ci{} defect.
Considering additionally the likewise lower diffusion barrier of \si{}, configurations of separated \cs{} and \si{} will occur in all probability.
-This is reinforced by the {\em ab initio} MD run at non-zero temperature, which shows structures of separating instead of recombining \cs{} and \si{} defetcs.
+This is reinforced by the {\em ab initio} MD run at non-zero temperature, which shows structures of separating instead of recombining \cs{} and \si{} defects.
This suggests increased participation of \cs{} already in the initial stages of the implantation process.
The highly mobile \si{} is assumed to constitute a vehicle for the rearrangement of other \cs{} atoms onto proper lattice sites, i.e. lattice sites of one of the the two fcc lattices composing the diamond structure.
-This way, stretched SiC strcutures, which are coherently aligned to the c-Si host, are formed by agglomeration of \cs.
+This way, stretched SiC structures, which are coherently aligned to the c-Si host, are formed by agglomeration of \cs.
Precipitation into an incoherent and partially strain-compensated SiC nucleus occurs once the increasing strain energy surpasses the interfacial energy of the incoherent 3C-SiC precipitate and the c-Si substrate.
As already assumed by Nejim~et~al.~\cite{nejim95}, \si{} serves as supply for subsequently implanted C atoms to form further SiC in the resulting free space due to the accompanied volume reduction.
%Indeed, combinations of \cs{} and \ci{} \hkl<1 0 0> DBs are observed.
%
Further conclusions are derived from results of the high C concentration simulations, in which a large amount of C atoms to obtain stoichiometry is incorporated into a small volume within a short period of time, which results in essentially no time for the system to rearrange.
-Due to this, the occurence of strong C-C bonds and the production of a vast amount of damage is observed, which finally results in the formation of an amorphous phase.
+Due to this, the occurrence of strong C-C bonds and the production of a vast amount of damage is observed, which finally results in the formation of an amorphous phase.
The strong bonds and damage obviously decelerate structural evolution.
The short time, which is not sufficient for structural evolution of the strongly damaged region, can be mapped to a system of low temperature, which lacks the kinetic energy required for the restructuring process.
This is assumed to be related to the problem of slow structural evolution encountered in the high C concentration simulations.
%
%Considering the efficiency of high implantation temperatures, experimental arguments exist, which point to the precipitation mechanism based on the agglomeration of \cs.
-More substantially, understoichiometric implantations at room temperature into preamorphized Si followed by a solid phase epitaxial regrowth step at \degc{700} result in Si$_{1-x}$C$_x$ layers in the diamond cubic phase with C residing on substitutional Si lattice sites \cite{strane93}.
+More substantially, understoichiometric implantations at room temperature into preamorphized Si followed by a solid-phase epitaxial regrowth step at \degc{700} result in Si$_{1-x}$C$_x$ layers in the diamond cubic phase with C residing on substitutional Si lattice sites \cite{strane93}.
The strained structure is found to be stable up to \degc{810}.
Coherent clustering followed by precipitation is suggested if these structures are annealed at higher temperatures.
%
% high t - epitaxial relation ...
Moreover, implantations below the optimum temperature for the IBS of SiC show regions of randomly oriented SiC crystallites whereas epitaxial crystallites are found for increased temperatures \cite{lindner99}.
The results of the MD simulations can be interpreted in terms of these experimental findings.
-The successive occupation of regular Si lattice sites by \cs{} atoms as observed in the high temperature MD simulations and assumed from results of the quantum-mechanical investigations perfectly statisfies the epitaxial relation of substrate and precipitate.
+The successive occupation of regular Si lattice sites by \cs{} atoms as observed in the high temperature MD simulations and assumed from results of the quantum-mechanical investigations perfectly satisfies the epitaxial relation of substrate and precipitate.
In contrast, there is no obvious reason for a topotactic transition of \ci{} \hkl<1 0 0> DB agglomerates, as observed in the low temperature MD simulations, into epitaxially aligned precipitates.
The latter transition would necessarily involve a much more profound change in structure.
% amorphous region for low temperatures
Indeed, the complex transformation of agglomerated \ci{} DBs as suggested by results of the low C concentration simulations could involve an intermediate amorphous phase probably accompanied by the loss of alignment with respect to the Si host matrix.
%
% perfectly explainable by Cs obvious hkl match but not for DBs
-In any case, the precipitation mechanism by accumulation of \cs{} obviously statisfies the experimental finding of identical \hkl(h k l) planes of substrate and precipitate.
+In any case, the precipitation mechanism by accumulation of \cs{} obviously satisfies the experimental finding of identical \hkl(h k l) planes of substrate and precipitate.
% no contradictions, something in interstitial lattice, projected potential ...
Finally, it is worth to point out that the precipitation mechanism based on \cs{} does not necessarily contradict to results of the HREM studies \cite{werner96,werner97,lindner99_2}, which propose precipitation by agglomeration of \ci.
\raisebox{600pt}{ }
\begin{tabular}{ll}
-Erstkorrektor: & Prof. Dr. Bernd Stritzker \\
-Zweitkorrektor: & Prof. Dr. Kai Nordlund \\
+Erstgutachter: & Prof. Dr. Bernd Stritzker \\
+Zweitgutachter: & Prof. Dr. Kai Nordlund \\
Tag der m"undlichen Pr"ufung: & 27. Juli 2011 \\
\end{tabular}
projects / lectures / latex.git / commitdiff