Base Doping Profile Control for SiGe PNP HBTs by Yiheng Lin B.Eng., Tsinghua University, 2008 M.A.Sc., The University of British Columbia, 2012 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in The Faculty of Graduate and Postdoctoral Studies (Materials Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) April 2016 © Yiheng Lin, 2016 ii Abstract The aim of this thesis is to investigate three aspects related to phosphorus diffusion for doping profile control in PNP SiGe HBTs: We systematically and quantitatively investigated the impact of carbon and Ge on P diffusion in strained SiGe:C up to 18% Ge and 0.32% C through experiments, which shows that the incorporation of carbon to retard P diffusion is not as effective in SiGe as it is in Si. Models were established to calculate the effective P diffusivities as a function of carbon concentration. These models can also be applied to boron, phosphorus, arsenic and antimony diffusion in Si with the presence of carbon. These results indicate that the microscopic mechanism of P diffusion in Si0.82Ge0.18 has a small but non-negligible vacancy-mediated term. An experimental study of thermal nitridation effects on phosphorus diffusion in strained Si1-xGex and strained Si1-xGex:Cy was performed. P diffusivities under thermal nitridation (vacancy injection) and the effective inert condition were compared. The result shows that thermal nitridation can retard P diffusion in SiGe with up to 18% Ge content, but the effectiveness of this retardation decreases with increasing Ge and C content. The Ge dependence can be explained by the increasing contribution from vacancy-assisted mechanism for P diffusion in strained SiGe with the increasing Ge content. P tends to segregate out of SiGe region, which happens simultaneously with diffusion. A coupled diffusion and segregation model is needed to predict the P profile evolution at thin SiGe layers. The model was re-derived theoretically, where the contributions from diffusion and segregation to dopant flux are explicitly shown. The model is generic to coupled diffusion and segregation in inhomogeneous alloys, and provides a new approach in segregation coefficient extraction. This model is especially helpful for heterostructures with lattice mismatch strains. Experiments of coupled P diffusion and segregation were performed with graded SiGe layers for Ge molar fractions up to 0.18, which are relevant to PNP SiGe HBTs. The model was shown to describe both diffusion and segregation behavior well. iii Preface Most of the work has been done in Department of Materials Engineering at UBC in the past three and half years, including the designs of experiments, thermal anneals with a Linkam heating stage, computer simulations with TSUPREM-4TM, Sentaurus ProcessTM and MATLABTM. After the wafer growth, sample preparation and patterning for various annealing conditions were performed by the author and was done mostly at Advanced Nanofabrication Facility (ANF) in Advanced Materials and Process Engineering Laboratory AMPEL at UBC. Experimentally, the author designed all the sample structures and annealing conditions, prepared samples for secondary ion mass spectrometer (SIMS). The X-ray diffraction (XRD) measurements were performed by the author using a PANalytical X’Pert PRO MRD tool at AMPEL. All the analysis of these characterizations, theoretical derivations, modeling and simulations were conducted by the author. Without the contributions from the collaborators, however, this thesis could not have been accomplished. All the samples were grown by Mr. Manfred Schiekofer at Texas Instrument Deutschland GmbH in Germany. The advanced rapid thermal anneals (RTAs) under inert nitrogen ambient were carried out by Mr. Bernhard Benna at Texas Instrument Deutschland GmbH in Germany. The advanced rapid thermal anneals (RTAs) under ammonia ambient were carried out by Mr. Ed Myers and Mr. Maurice Stevens at Stanford Nanofabrication Facility at Stanford University. The SIMS measurements were performed by Evans Analytical Group (EAG), a provider of commercial material characterizations. In addition, two journal articles have been published and another journal article manuscript has been submitted based on the work in Chapter 4 to 6. 1. Y. Lin, H. Yasuda, M. Schiekofer, B. Benna, R. Wise, and G. (Maggie) Xia, “Effects of carbon on phosphorus diffusion in SiGe:C and the implications on phosphorus diffusion mechanisms,” Journal of Applied Physics, vol. 116, no. 14, p. 144904, 2014. iv 2. Y. Lin, H. Yasuda, M. Schiekofer, and G. (Maggie) Xia, “Coupled dopant diffusion and segregation in inhomogeneous SiGe alloys: Experiments and modeling,” Journal of Applied Physics, vol. 117, no. 21, p. 214901, 2015. 3. Y. Lin, H. Yasuda, M. Schiekofer, and G. (Maggie) Xia, “The Effect of Thermal Nitridation-Based Defect Engineering on Phosphorus Diffusion in SiGe:C,” Journal of Materials Science, 51(3), 1532-1540, 2015. v Table of Contents Abstract ........................................................................................................................................... ii Preface ............................................................................................................................................ iii Table of Contents ............................................................................................................................ v List of Tables ................................................................................................................................ viii List of Figures ................................................................................................................................ ix List of Abbreviations .................................................................................................................... xiii Acknowledgements ...................................................................................................................... xiv Dedication .................................................................................................................................... xvi Introduction ..................................................................................................................... 1 1.1 Industry background, motivations and problem definitions ................................................. 1 1.1.1 Heterojunction bipolar transistor basics ........................................................................ 1 1.1.2 Requirement of narrow base doping profiles for HBTs ................................................ 3 1.2 Problem definitions ............................................................................................................... 5 1.3 The organization of this thesis .............................................................................................. 7 Fundamentals of dopant diffusion in semiconductors .................................................... 9 2.1 Macroscopic description of diffusion ................................................................................. 10 2.2 Microscopic description of diffusion .................................................................................. 11 2.2.1 Atomic-scale description of diffusion in solids ........................................................... 11 2.2.2 Point defects and their roles in dopant diffusion ......................................................... 12 2.2.2 Fermi level effect in dopant diffusion ......................................................................... 15 2.3 Methods of dopant profile engineering .............................................................................. 16 Experiment and modeling methodologies .................................................................... 18 3.1 General procedure of diffusion study in this work ............................................................. 18 3.2 Sample structure design and growth ................................................................................... 19 3.2.1 Strain and critical thickness for SiGe epitaxial heterostructures ................................. 19 3.2.2 Sample structure and growth ....................................................................................... 21 3.3 Annealing of diffusion samples .......................................................................................... 21 3.4 Dopant profile characterization: secondary ion mass spectrometry ................................... 22 3.5 Strain characterization by X-ray diffraction ....................................................................... 23 3.5.1 Principles of X-ray diffraction ..................................................................................... 23 vi 3.5.2 Reciprocal space mapping ........................................................................................... 25 3.6 Computer modeling and simulations .................................................................................. 26 3.6.1 Simulations with process simulator TSUPREM-4TM .................................................. 27 3.6.2 Finite difference method .............................................................................................. 28 The effect of carbon on phosphorus diffusion in SiGe ................................................. 31 4.1 Introduction ........................................................................................................................ 31 4.2 Literature review: the effect of carbon on dopant diffusion in Si and SiGe ....................... 32 4.2.1 Carbon in Si and SiGe ................................................................................................. 32 4.2.2 Substitutional carbon incorporation in epitaxial grown Si1-yCy and Si1-x-yGexCy......... 33 4.2.3 Mechanisms of carbon in defect and diffusion engineering ........................................ 36 4.2.4 Carbon impacts B and P diffusion in Si similarly ....................................................... 37 4.2.5 Carbon impacts B and P diffusion in SiGe differently ................................................ 38 4.3 Diffusion experiments and results ...................................................................................... 38 4.3.1 Sample structures and annealing conditions design .................................................... 38 4.3.2 Strain and depth profile analysis ................................................................................. 40 4.3.3 P diffusivity ratio extraction ........................................................................................ 44 4.4 Modeling of carbon effect on P diffusion ........................................................................... 47 4.5 Model parameter extraction ................................................................................................ 52 4.6 Comparison of carbon impact on different dopants ........................................................... 53 4.7 Conclusions ........................................................................................................................ 56 The effect of thermal nitridation on phosphorus diffusion in SiGe and SiGe:C ........... 57 5.1 Literature review ................................................................................................................ 57 5.1.1 Physical effects of thermal nitridation ......................................................................... 57 5.1.2 Effect of thermal nitridation on dopant diffusion ........................................................ 59 5.1.3 Thermal nitridation on Si and SiGe ............................................................................. 61 5.2 Diffusion experiments ........................................................................................................ 61 5.2.1 Sample structures and annealing conditions design .................................................... 62 5.2.2 Defect injection using masking layers ......................................................................... 63 5.3 Experimental results ........................................................................................................... 64 5.3.1 Strain and composition analysis .................................................................................. 64 5.3.2 Concentration profiling and P diffusivity ratio extraction ........................................... 66 5.4 Discussions ......................................................................................................................... 72 5.4.1 Major observations ...................................................................................................... 72 vii 5.4.2 Ge fraction dependence of thermal nitridation impact factor ...................................... 72 5.4.3 Carbon fraction dependence of thermal nitridation impact factor ............................... 73 5.4.4 Thermal nitridation and equivalent [𝑰] 𝑰 ∗ in Si0.82Ge0.18 ............................................ 74 5.5 Conclusions ........................................................................................................................ 77 Segregation of phosphorus in SiGe............................................................................... 78 6.1 Introduction ........................................................................................................................ 78 6.2 Literature review ................................................................................................................ 80 6.3 Derivation of coupled of diffusion-segregation equation ................................................... 83 6.3.1 The driving force for diffusion in solid state materials ............................................... 84 6.3.2 Derivation of coupled diffusion-segregation equation ................................................ 84 6.4 Experimental verification ................................................................................................... 88 6.4.1 Sample structure design and experiments ................................................................... 88 6.4.2 Parameters for modeling .............................................................................................. 89 6.4.3 SIMS data and data fitting ........................................................................................... 91 6.4.4 Model prediction and P profile time evolution ............................................................ 95 6.5 Impact of P segregation on HBT device performance ....................................................... 98 6.6 Conclusions ...................................................................................................................... 101 Contributions and future work .................................................................................... 103 7.1 Contributions of this work ................................................................................................ 103 7.2 Suggestions for future work ............................................................................................. 104 Bibliography ................................................................................................................................ 106 Appendix A Numerical fitting with TCAD softwares ................................................................ 112 Appendix A1 P diffusivity ratio extraction in TSUPREM-4 ................................................. 112 Appendix A2 Carbon/Ge model implementation in TSUPREM-4 ........................................ 116 Appendix A3 Carbon/Ge model implementation in Sentaurus Process ................................. 121 Appendix B Numerical Fitting with MATLAB ......................................................................... 125 Appendix B1 Numerical solver for diffusion-segregation equation....................................... 125 Appendix B2 Numerical search for best fitting with least-square method ............................. 128 Appendix C Remarks from the author ........................................................................................ 130 viii List of Tables Table 2.1 Terminology used for diffusion conditions. ..............................................................................14 Table 5.1 Epitaxial wafer matrix used in this work. ................................................................................62 Table 5.2 Nitridation impact factor for phosphorus RNTD,P under non-equilibrium point defect conditions caused by nitridation. The annealing temperatures quoted here are the temperature extracted using TSUPREM-4 fitting to account for uncertainty. ...................................................71 Table 5.3 Summary of the impact factor of carbon, thermal nitridation, and both methods on P diffusion in Si0.82Ge0.18 and the corresponding interstitial undersaturation ratio [I]/[I]*. The impact of carbon is adopted from Ref. [91]. .....................................................................................75 Table 6.1 List of major parameters used in TSUPREM-4 simulations. The diffusivities are in μm2/min. Concentrations are in cm-3. ................................................................................................................90 ix List of Figures Figure 1.1 Opposite doping configurations in NPN and PNP type BJTs. Figure is reprinted from Wikipedia under the GNU Free Documentation License. ................................................................ 2 Figure 1.2 NPN BJT with forward-biased E–B junction and reverse-biased B–C junction. Figure is reprinted from Wikipedia under the GNU Free Documentation License. ..................................... 2 Figure 1.3 Cross-section view of a SiGe HBT [1]. ...................................................................................... 3 Figure 1.4 Dopant and Ge profiles in the design of 800/1000 GHz fT/fmax SiGe HBTs [3-5]. ................. 4 Figure 2.1 Illustration of different types of point defects and dopant-defect pairs in crystal. ..............13 Figure 3.1 The standard procedure of this work. .....................................................................................18 Figure 3.2 The critical thickness vs. Ge fraction for pseudomorphic Si1−xGex layers grown on bulk (100) silicon . Figure reprinted from Ref.[1] with permission. .................................................................19 Figure 3.3 A cross-sectional illustration of Bragg condition. ..................................................................24 Figure 3.4 Two dimensional reciprocal lattice and the Ewald sphere representation. In a reciprocal space map, the alignment of the diffraction spot from a SiGe layer with respect to the spot from the underlying Si substrate shows whether this layer is fully strained (i.e., with the same in-plane lattice constant) with the Si substrate. ..............................................................................................25 Figure 3.5 Illustration of finite difference method for solving mass transport equations. ...................28 Figure 3.6 Illustration showing how diffusion flux such as Fick’s first law is expressed by finite difference equations. ...........................................................................................................................29 Figure 3.7 Illustrations of numerical calculation procedures across spatial and time domain............30 Figure 4.1 Summary of the carbon substitutionality in epitaxial Si1-yCy for different growth temperatures and silane partial pressures. Figure reprinted from Ref.[60] with permission. .....35 Figure 4.2 Summary of carbon substitutionality in epitaxial Si1-x-yGexCy for different growth temperatures and carbon sources. The Si1-x-yGexCy samples include Ge contents from 10 to 36%. Figure reprinted from Ref.[63] with permission. .............................................................................36 Figure 4.3 Schematic diagram of the structure design and splits in C mole fractions. The dotted line indicates the original wafer surface. The carbon concentration split is from 0.06 to 0.32 at.%, x which corresponds to 3e19 to 1.6e20 cm-3. The Ge and C profiles are designed to be constant with a certain thickness as shown above. ..................................................................................................40 Figure 4.4 (113) reciprocal space mapping of annealed samples with the largest thermal budget for each type of structures. a) Ge18C000; b) Ge18C009; c) Ge18C018; d) Ge18C032. ......................41 Figure 4.5 XRD (004) scan of SiGe and SiGe:C as-grown samples showing carbon molar fraction change. .................................................................................................................................................43 Figure 4.6 As-grown and annealed Ge, P and C SIMS profiles of Ge18C018 samples. ........................44 Figure 4.7 Example of TSUPREM-4 fitting with SIMS profiles of as-grown and annealed samples Ge18C018. ...........................................................................................................................................46 Figure 4.8 Summary of the carbon impact factor for P and B diffusivity in Si1-xGex:Cy. Literature data are from Ref.[12], [15], [16], [48]. ......................................................................................................46 Figure 4.9 The relation between the interstitial suppression ratio ([I]/[I]*) and dopant diffusivity suppression ratio (RC = D/D*) for different diffusion mechanism under the assumption of [I][V] = [I]*[V]*. .................................................................................................................................................50 Figure 4.10 Simplified illustration of cross-section view of Si1-xGex and Si1-xGex:Cy lattice spacing under different strain status: (a) Si and fully strained Si1-xGex:Cy with x/y ≈ 10, and the out-of-place lattice constant is the same as the Si out-of-plane lattice constant, (b) fully strained Si1-xGex:Cy with x/y > 10 , (c) fully strained Si1-xGex and (d) fully relaxed Si1-xGex . .......................................51 Figure 4.11 Modeling of carbon effect on P and B diffusivity in Si and strained SiGe. Literature data are from the same references as in Figure 4.8. .................................................................................53 Figure 4.12 Experiment RC data for Sb, As, P and B vs. model predictions. Literature data are from the same references as in Figure 4.8. .................................................................................................55 Figure 5.1 An example of the etch figures of the stacking faults taken in a wafer after thermal nitridation at 1100 ℃ for 4 hours [86]. .............................................................................................58 Figure 5.2 The shrinkage of stacking faults during nitridation of silicon without oxide film under an ammonia pressure of 1 kg/cm2 [86]. .................................................................................................59 Figure 5.3 Schematic diagram of the structure design. ...........................................................................63 xi Figure 5.4 Structures used for different defect injection conditions during thermal nitridation experiments. ........................................................................................................................................64 Figure 5.5 (113) reciprocal space mapping of annealed samples with the largest thermal budget for each type of structures: a) Ge09C0; b) Ge09C009; c) Ge18C0; d) Ge18C009. .............................65 Figure 5.6 XRD (004) scan of SiGe and SiGe:C samples used for thermal nitridation experiments. The shifts of (004) peaks show Ge and carbon molar fraction differences. ...........................................66 Figure 5.7 SIMS profiles of as-grown and annealed samples and TSUPREM-4 fittings. The solid lines on top of the P profiles are the TSUPREM-4 fittings. .....................................................................69 Figure 5.8 Nitridation impact factor RNTD as a function of carbon % at various Ge %. Error bars of the data without carbon are not shown for clarity. The literature data quoted in this figure are obtained from Ref.[36] for P in Si and from Ref.[16] for B in Si and SiGe. ..................................71 Figure 5.9 P diffusivity under a defect engineering condition divided by the P diffusivity ratio as a function of the interstitial undersaturation ratio [I]/[I]* in Si0.82Ge0.18. It is based on local equilibrium assumption of [I][V] = [I]*[V]*.This figure was adapted from Ref. [91]. ...................76 Figure 6.1 Illustration of sharp and graded interfaces. ...........................................................................79 Figure 6.2 Comparison of phosphorus concentration profiles in Ge-doped (solid curve) and control (dotted curve) silicon. Figure reprinted from Ref. [98] with permission. ......................................81 Figure 6.3 Schematic diagram of the structure design. The Ge peak concentrations used were 14% and 18%. .....................................................................................................................................................88 Figure 6.4 Experimental data of P segregation coefficient in SiGe. The data points in diamonds are from Ref. [92], which were obtained from SiGe films thicker than 0.2 micron where no strain analysis was available. The remaining data are extracted from Ref. [15], [74], [98]. ...................91 Figure 6.5 Ge and P SIMS profiles and best-fitting curves using the coupled diffusion-segregation model in Eq. (6.18). Simulations with diffusion-only model were plotted for comparison. ..........94 Figure 6.6 The best-fitting Eseg and the corresponding segregation coefficients extracted from coupled diffusion-segregation experiments. The dashed parts are extrapolation using the Eseg values. Literature data are from the same reference as in Figure 6.4. .......................................................95 xii Figure 6.7 (a) Simulated time evolution of P profiles during annealing in comparison with SIMS data. The peak Ge fraction is 14%, (b) simulated time evolution of P profiles across sharp Si/SiGe interfaces during annealing and (c) extracted P concentration ratio at either side of the left Si/SiGe interface as a function of time (from 1 to 10 hours). The initial P and Ge have box-like profiles. The kseg for Si/Si0.8Ge0.2 (fully strained) is shown for comparison. ...................................98 Figure 6.8 Doping and Ge profiles used for one-dimensional HBT device simulation. ......................100 Figure 6.9 One-dimensional device simulation results showing the impact of the P segregation effect on PNP HBT current gain (𝜷). The P base profiles are simulated using the diffusion-segregation model and the diffusion only model for comparison. The annealing conditions used are 5 s and 10 s at 1000 ⁰C. .......................................................................................................................................100 Figure 6.10 One-dimensional device simulation results showing the impact of the P segregation effect on PNP HBT cut-off frequency (fT). The P base profiles are simulated using the diffusion-segregation model and the diffusion only model for comparison. The annealing conditions used are 5 s and 10 s at 1000 ⁰C................................................................................................................101 xiii List of Abbreviations ALD Atomic Layer Deposition BJT Bipolar Junction Transistor CVD Chemical Vapor Deposition CMOS Complementary Metal-Oxide-Semiconductor C-BiCMOS Complementary BJT-CMOS FEA Finite Element Analysis HBT Heterojunction Bipolar Transistor HDD Hard Disk Drive ITRS International Technology Roadmap for Semiconductors ODD Optical Disk Drive RTA Rapid Thermal Annealing RTP Rapid Thermal Processing RSM Reciprocal Space Mapping SiGe Silicon-Germanium SIMS Secondary Ion Mass Spectrometry TFT Thin Film Transistor TCAD Technology Computer Aided Design XRD X-ray Diffraction xiv Acknowledgements First and foremost, I would like to thank my thesis advisor, Prof. Guangrui (Maggie) Xia, for her guidance, patience and trust throughout my graduate journey. Without her enduring support and encouragement, it would not have been possible for me to keep on pursuing a doctorate degree. She not only taught me how to conduct high-quality research, but also shared with me many experiences on how to be a professional engineer and work with the industry. I would also like to thank Prof. David Dixon and Prof. Chad Sinclair for being on my thesis supervisory committee and Prof. Warren Poole for being the chairman of my oral examining committee. Their insightful questions and constructive advice helped me a lot in improving my work. I am very grateful to Texas Instruments for funding this project through Semiconductor Research Corporation. I would like to thank our collaborators at Texas Instruments (TI): Hiroshi Yasuda for his constant support and kind help during the entire project, Dr. Stanley Philips for helpful discussions and advice, Mr. Manfred Schiekofer and Mr. Bernhard Benna at TI Germany for growing all the wafers and performing the annealing experiments. I would also like to thank Dr. Rick Wise and Dr. Howard Ho for facilitating the funding process from TI. I would like to thank Dr. Mario Beaudoin and Dr. Alina Kulpa at Advanced Nanofabrication Facility at UBC for their training and generous help in many aspects of my cleanroom experience. I would like to thank Dr. Maurice Steven and Dr. Ed Myers at Stanford Nanofabrication Facility for their help on annealing experiments. I would like to thank Dr. Bob Hengstebeck at Evans Analytical Group for helpful discussions and advice on SIMS analysis, which is crucial to this work. I am indebted to Prof. James Plummer at Stanford University and Prof. Scott Dunham at University of Washington for helpful discussions and advice. xv I would like to give my special thanks to Ms. Michelle Tierney, Ms. Mary Jansepar, and Ms. Fiona Webster at UBC Materials Engineering for their kind help on many administrative issues. I would like to thank Dr. Yuanwei Dong, my former lab mate, for many brain-storming discussions and for proofreading my first paper. I would also like to thank my current lab mates Mr. Feiyang Cai, Ms. Ye Zhu, Mr. Xiyue Li and Mr. Weijun Luo for their kind help. I would like to thank Mr. Aron Fredrick for proofreading of my thesis. Finally, I wish to thank my family, my parents, my wife and my son, for their ever-lasting love and faith on me. “Not everything that counts can be counted.” I always owe them for their unconditional love. xvi Dedication To my parents, who gave and nurtured my life. 献给我的父母,是他们给予并哺育我的生命。 To my wife, who tolerates, inspires, enriches, and passes on my life. 献给我的妻子,是她容忍、启迪、丰富并传递我的生命。 To my son, who will carry on my life into the future. 献给我的儿子,是他将延续我的生命直到未来。 To nature, the origin and destination of my life. 献给大自然,我的原点和归宿。 1 Introduction 1.1 Industry background, motivations and problem definitions For decades, the silicon (Si) based semiconductor industry has been crucial to drawing the broad outlines of the modern history of economics and technology. Silicon has been the dominant material of the semiconductor industry for over 40 years. Thanks to its many excellent material properties, low cost, and superior manufacturability, silicon based products account for over 97% of the microelectronics market [1]. As close relatives to Si, Silicon-Germanium (SiGe) and SiGe:C alloys emerged and have been applied widely in the microelectronics industry for the past two to three decades. This is due to their capability in strain and band gap engineering and the compatibility to the mainstream Si processing. The high speed SiGe heterojunction bipolar transistor (HBT) is one such example. Since the first sale in 1998, SiGe HBTs have seen increasing applications in wireless communications, hard disk drive (HDD), optical disc drive (ODD) for CD, DVD and Blu-ray Disk, thin-film transistor (TFT) displays, high resolution video, optical networking components and so on. 1.1.1 Heterojunction bipolar transistor basics A heterojunction bipolar transistor (HBT) is a type of bipolar junction transistor (BJT). BJTs are made of two p-n junctions and come in two types, or polarities, known as PNP and NPN. The type is determined by the doping of the three main terminal regions: emitter, base and collector (Figure 1.1). Bipolar transistors are so named because their operation involves both electrons and holes. Charge flow in a BJT is due to bi-directional flows of charge carriers across a junction between two regions of different charge concentrations. Figure 1.2 shows the normal operation of a BJT where the emitter-base (E–B) junction is forward-biased and the 2 base-collector (B-C) junction is reverse-biased. It also shows how the carriers transport during operation. Figure 1.1 Opposite doping configurations in NPN and PNP type BJTs. Figure is reprinted from Wikipedia under the GNU Free Documentation License. Figure 1.2 NPN BJT with forward-biased E–B junction and reverse-biased B–C junction. Figure is reprinted from Wikipedia under the GNU Free Documentation License. For high-speed and high-frequency applications, there are two important figures-of-merit: the common-emitter short-circuit cut-off frequency 𝑓𝑇 and the maximum frequency of oscillation 𝑓𝑚𝑎𝑥. Higher 𝑓𝑇 and 𝑓𝑚𝑎𝑥 represent better HBT device performance. 3 Traditionally, the emitter, the base and the collector of a BJT are fabricated in the same bulk material, such as Si. If different materials are used for emitter and base regions, a heterojunction is created at the emitter-base junction and this BJT becomes an HBT. Compared to homojunction BJTs, HBTs have several advantages: such as higher 𝑓𝑇 and 𝑓𝑚𝑎𝑥, improved injection efficiency and higher current gain, etc. Such characteristics are highly desired in many modern systems, especially for high-speed analog and mixed-signal (digit-analog) circuits. HBTs were dominated by III-V compound semiconductors before the emergence of SiGe HBTs, which were enabled by the lattice-matched strained-SiGe epitaxy technology in the 1990s. Figure 1.3 shows an example of a modern NPN SiGe HBT structure, where a SiGe layer is grown epitaxially on a Si substrate [2]. Figure 1.3 Cross-section view of a SiGe HBT [1]. 1.1.2 Requirement of narrow base doping profiles for HBTs Reducing the base width is an effective way to reduce the base transit time and thus increase the 𝑓𝑇. Shi et al. performed two-dimensional device simulations for the design of 800/1000 GHz 𝑓𝑇/𝑓𝑚𝑎𝑥 SiGe HBTs and shows that very aggressively scaled (high concentration and narrow width) vertical doping profile in the base region is required to achieve this high frequency [3], [4]. Physically, higher base doping reduces electrical resistance, while narrower base doping profile Horizontal Vertical 4 reduces capacitance and carrier transit time in the base region. These all contribute to a higher 𝑓𝑇/𝑓𝑚𝑎𝑥 and better high-frequency performance. Bolze et al. recently reported experimental results of 𝑓𝑇/𝑓𝑚𝑎𝑥 for up to 220/280 GHz which resulted from an aggressively scaled base doping profile and an advanced annealing technique, which is shown in Figure 1.4 [5]. Although the above simulations and experimental demonstrations were performed on NPN SiGe HBTs, the performance of PNP type HBTs benefit from a well-controlled based doping profile in the same way. Figure 1.4 Dopant and Ge profiles in the design of 800/1000 GHz fT/fmax SiGe HBTs [3-5]. Compared to NPN SiGe HBTs, PNP SiGe HBTs are much less explored due to the inherently slower minority carrier transport and thus lower frequency response. However, complementary SiGe HBTs (with both NPN and PNP HBTs) have many advantages over an NPN-only technology for numerous analog applications requiring high speed, low noise, and large Ge Doping 5 voltage swing [6], [7]. Furthermore, complementary-SiGe bipolar complementary metal-oxide-semiconductor (SiGe C-BiCMOS) technology platform containing both NPN and PNP devices with, ideally, equal performance, and CMOS can offer compelling advantages in many types of analog circuits. It can be used for lowering voltage supply rails, improving current sources, and implementing fast and low distortion driver amplifiers. To date, more and more effort is involved in the R&D of complementary-SiGe BiCMOS technology (C-BiCMOS), which utilize both NPN and PNP HBTs [8], [9]. The International Technology Roadmap for Semiconductors (ITRS) 2013 edition on radio frequency (RF) and Analog/Mixed-signal technologies for wireless communications has addressed the needs for high speed PNP HBTs (HS-PNP), which is driven by applications that require high performance analog and mixed-signal ICs [10]. These ICs, such as drivers, video amplifiers, bus interfaces, digital/analog and analog/digital converters (DAC/ADC), and operational amplifiers greatly benefit from C-BiCMOS technologies that offer both NPN and PNP transistors with matched performances (comparable frequency performance). Matching the performance of the NPN transistor is the aim of the HS-PNP roadmap. Based on the above information, we can see that aggressively scaled vertical doping profiles are essential in improving HS-PNP HBTs’ performance. The state-of-the-art PNP HBTs for C-BiCMOS technologies use phosphorus (P) as the base layer dopant and strained SiGe as the base layer material. Therefore, to achieve narrow base doping, it is a must to understand and predict P diffusion and segregation in strained SiGe, which is the goal of this thesis work. 1.2 Problem definitions To achieve this goal, we need to first introduce major approaches to suppress or limit phosphorus (P) diffusion in Si and SiGe. To date, there are two major approaches for this purpose: 1. Defect engineering methods such as carbon incorporation or thermal nitridation. 6 2. Advanced thermal processing techniques such as sub-millisecond annealing techniques with less thermal budget and therefore less diffusion. This work is not about option 2, as rapid thermal processing (RTP) is the industry mainstream thermal processing technique. Existing studies of P diffusion in strained-SiGe have been performed mostly using furnace annealing, mainly because of the availability of the tools, long anneal time and better temperature accuracy for diffusion study. Studies using RTP tools are in great need to be industry relevant. About option 1, studies have shown that carbon incorporation in Si and strained-SiGe can cause the unsaturation of Si self-interstitials [11], [12], therefore interstitial-mediated dopants are suppressed, such as boron [13], [14] and phosphorus diffusion [13], [15]. Karunaratne et al. reported the first quantitative comparison of B diffusivities in Si, Si0.9Ge0.1, Si0.999:C0.001 and Si0.899Ge0.1:C0.001 to extract carbon retardation factors [16]. It has been a common industry practice to use carbon in SiGe to retard B diffusion for NPN SiGe HBTs while there are only a handful of studies on P diffusion in SiGe:C for PNP HBTs. Tillack et al. studied carbon retardation effect in Si0.8Ge0.2 on P diffusion for two different carbon concentrations after annealing at 850 ⁰C [9]. However, this temperature is far from the temperature range where most diffusion happens in industry practice, which is above 1000 ⁰C. Systematic and quantitative study of carbon retardation effect on P diffusion in SiGe is still lacking. This is highly desired in industry for accurate device and process design of HS-PNP HBTs and will be addressed as Topic 1 of this study, and will be discussed in Chapter 4. The results in topic 1 are also needed to serve as the base line in the study of thermal nitridation effect in Topic 2. Thermal nitridation is also known to suppress B and P diffusion in Si via exposing bare Si surface in ammonia ambient, where vacancies are injected into Si substrates and retard interstitial diffuser’s motion [17]. This effect was studied in late 80’s and early 90’s mainly for diffusion mechanism study [18], [19], and has no published record on its industry applications. R&D on SiGe 7 PNP-HBTs provides a good opportunity for the application of this method, as the dopants in PNP-HBTs, B and P, are both interstitial diffusers, and can be retarded at the same time with vacancy injection method. However, there is no study available on whether this method also works in strained-SiGe system and (if any) to what extent. It is also not clear whether this effect will work with other approach (e.g. carbon incorporation) and to what extent. If this method is effective, it is a relatively low cost and simple approach to be adopted by industry, as the only change is in the annealing ambient. These questions will be address as Topic 2 in Chapter 5. Topic 3 is about P segregation across graded SiGe layers (common SiGe layer design for HBTs), which happens simultaneously with P diffusion. As P prefers to stay in Si rather than in SiGe, this effect is not desired for narrow P profiles, and needs to be taken into account in the mass transport modeling. While several studies on this effect are available, but none of them are for graded SiGe layers and lack experimental verification and appropriate modeling. This topic will be discussed in Chapter 6. 1.3 The organization of this thesis After this chapter of introduction, Chapter 2 sets the stage for the remainder of this thesis by presenting the fundamental concepts of dopant diffusion in semiconductors. In particular, we discuss the relationship between macroscopic impurity diffusion and the underlying microscopic diffusion mechanisms in elemental semiconductors like silicon. Chapter 3 presents the methodologies and procedures of the experiments and computer modeling used in this work. Then, three topics of our study related to the P doping profile control are presented in Chapters 4 to 6. In Chapter 4, we present the experiments and modeling of carbon impact on P diffusion in strained SiGe for up to18% Ge and 0.32% carbon. The impact of carbon and Ge on P diffusion in strained SiGe:C was systematically and quantitatively investigated by experiments. Models were established to calculate the effective P diffusivities as a function of carbon concentration, which 8 can also be applied to boron, phosphorus, arsenic and antimony diffusion in Si with the presence of carbon. These results indicate that the microscopic mechanism of P diffusion in strained SiGe may have a small but non-negligible vacancy-mediated term. In Chapter 5, we present the experiments of thermal nitridation on P diffusion in strained SiGe. The possibility of applying thermal nitridation as a technology to suppress P diffusion in strained SiGe and SiGe:C is explored. In Chapter 6, we present the experiments and modeling of P segregation across graded SiGe layers. A coupled diffusion and segregation model is re-derived based on thermodynamic principles, where the contributions from diffusion and segregation to dopant flux are explicitly shown. Experiments of coupled P diffusion and segregation were performed with graded SiGe layers for Ge molar fractions up to 0.18, which are relevant to PNP SiGe HBTs. It is shown that this model can well describe the diffusion and segregation behavior simultaneously. The model is generic to coupled diffusion and segregation in inhomogeneous alloys. It also provides a new approach in segregation coefficient extraction, which is especially helpful for heterostructures with lattice mismatch strain. Chapter 7 summarizes the contributions of our work and suggests relevant topics to further explore in the future. 9 Fundamentals of dopant diffusion in semiconductors This thesis is concerned with a technological ambition to find ways to slow down P diffusion in strained SiGe and to predict P diffusion profiles. This task requires the understanding of diffusion and the use of various phenomenological and atomistic theories. On a macroscopic scale, diffusion can, in most cases, be described by Fick’s law, a simple second order partial differential equation. On a microscopic scale, however, it is a complicated many-body physics problem, whose exact solutions are almost impossible to calculate. In between these two extremes, there is a whole area where most of today’s diffusion research is focused on. First-principles and atomistic calculations attempts to work from the microscopic end up, while diffusion experiments along with continuum simulations, such as the work in this thesis, try to work from the macroscopic end down. Unfortunately, the gap between the macroscopic and microscopic understanding of solid-state diffusion is still waiting to be converged [20]–[23]. In other words, there is not a one-to-one correspondence between a microscopic diffusion mechanism and a macroscopic impurity diffusivity. The goal of this chapter is to introduce them as a basis for more in-depth discussions in the following chapters. When using physical or scientific concepts to solve real-world problems, it is always beneficial to keep in mind that “any serious consideration of a physical theory must take into account the distinction between the objective reality, which is independent of any theory, and the physical concepts with which the theory operates. These concepts are intended to correspond with the objective reality, and by means of these concepts we picture this reality to ourselves” [24]. In fact, if we parallel the history of atomic theory and diffusion studies, one can realize that the research on diffusion (since early 19th century) proceeded ahead of the acceptance of atomic theory (at around early 20th century), without a clear knowledge of what was diffusing exactly [25]. 10 2.1 Macroscopic description of diffusion Diffusion is the transport of matter from one point to another. The Latin word “diffundere” means “to spread out” [26]. Depositing a droplet of ink in a basin of water without stirring gives a simple demonstration of diffusion. After a period of time, the solution will become uniformly coloured. Adolf Fick proposed a mathematical framework to describe diffusion phenomena macroscopically [27]. It was postulated that the flux of matter 𝑗 in the 𝑥 direction is proportional to the pertaining gradient of concentration 𝐶: 𝑗 = −𝐷𝜕𝐶𝜕𝑥 (2.1) This is what we call Fick’s first law today. 𝐷 is denoted as the diffusion coefficient, or diffusivity. In his paper Fick described it as “a constant dependent upon the nature of the substance”, which depicts that how fast a substance moves. Using the principle of conservation of mass, Fick also derived a second law of diffusion 𝜕𝐶𝜕𝑡=𝜕𝜕𝑥(𝐷𝜕𝐶𝜕𝑥) (2.2) If we look at Fick’s original paper, one may notice that the above equations have been updated and optimized compared to his original expression 𝛿𝑦𝛿𝑡= 𝑘𝛿2𝑦𝛿𝑥2 (2.3) where 𝑦 is the concentration (an extra minus sign is ignored here). There was an implicit assumption behind this equation that diffusivity is independent of concentration and of its gradient. When it comes to solid state matters, the atoms are closely stacked with certain pattern (for example, crystal structure) and their mutual interactions impact their motion. The diffusivity of a diffusing species in solid state matters represents not only the nature of the species itself, but also its interaction with the diffusion media. Fick’s provide only a phenomenological approach and respective empirical mathematical solutions from a macroscopic perspective. For a microscopic viewpoint, however, we would like 11 to understand how the individual atoms (the “substance” mentioned by Fick) move from one point to another and how their behavior finally sums up towards the macroscopic phenomenon. No matter whether we would like to harness (e.g. the drive-in diffusion for the introduction of dopants) or prevent (e.g. the out-diffusion of dopants during unavoidable high temperature processes, as studied in this work) diffusion, such understandings can help to find suitable technological solutions. 2.2 Microscopic description of diffusion In a crystalline solid, the crystal lattice restricts the positions and the migration paths of a diffusing species. This feature contrasts with a gas, where random distribution and displacements of atoms are assumed, and with liquids and amorphous solids, which are neither really random nor really ordered. The obstacle for the investigation of diffusion in solids is that it is extremely challenging to trace the motion of individual atoms in real time, even with modern high-resolution electron microscopy. Therefore the understanding of atomistic diffusion mechanism comes from a combined interpretation of macroscopic experimental results and microscopic diffusion theory. 2.2.1 Atomic-scale description of diffusion in solids Above zero Kelvin, every atom or ion in a solid material is vibrating very rapidly about its lattice position within the crystal at a frequency on the order of 1013 𝐻𝑧 with a typical vibration amplitude of about 10−12𝑚 [28]. For simplicity, we will use “atoms” to stand for the elementary particles that constitute a solid diffusion media material instead of using “atoms or ions”, as here the diffusion media are Si and SiGe. In a piece of material, there are numerous atoms whose kinetic energies are not equal and follow a natural and dynamic distribution. Thus, there are always probabilities for the atoms with higher energy to be too active to stay in their original positions. Given some time, these atoms will 12 relocate themselves and find new positions. In other words, it is this “hopping game” which finally constitutes the diffusion towards a statistically macroscopic scale. In a nearly perfect crystal material, the most simple atomic diffusion mechanism for substitutional impurities would be to exchange places with adjacent atoms [23]. Such a possibility has been investigated by Huntington and Seitz in their pioneering studies of two-atoms direct exchange in copper [19, 20]. The results indicate that a direct interchange is not favored in terms of energetic considerations. This diffusion mechanism which does not involve any defects has been investigated all the way until today, mainly from a first-principles perspective [31]. However, both theoretical and experimental results suggest that the contribution from exchange mechanism can be ignored [21–23]. On the other hand, the assistance of point defects in diffusion is believed to be of dominant importance in most crystalline materials, such as metals and semiconductors [24, 25]. 2.2.2 Point defects and their roles in dopant diffusion In crystalline materials, there are many types of defects ranging from zero to three dimensional. For crystals with low defect densities such as crystalline semiconductors, it is the point defects (zero-order) which act as the vehicles for diffusion. Figure 2.1 shows the major types of point defects. In Si crystals, for example, a vacancy is a missing atom from lattice and self-interstitial is an extra Si atom in positions other than lattice site. These two types of point defects are called native point defects, or intrinsic point defects. Impurity atoms are another type of point defects, which can occupy lattice sites or interstitial positions in the lattice. In the following discussions, we refer to native point defects simply as point defects unless otherwise noted. 13 Figure 2.1 Illustration of different types of point defects and dopant-defect pairs in crystal. The idea that dopant diffusion in silicon happens through interstitial-assisted and vacancy-assisted mechanisms has now been generally accepted based on both experimental observations and theoretical calculations [17, 26]. This concept can be formulated by 𝐷𝐴𝑒𝑓𝑓 = 𝐷𝐴∗(𝑓𝐼[𝐼] [𝐼]∗+ 𝑓𝑉[𝑉] [𝑉]∗) (2.4) where 𝐷𝐴𝑒𝑓𝑓 is the effective diffusivity of the dopant A measured under conditions where the point-defect populations are perturbed and 𝐷𝐴∗ is the normal, equilibrium diffusivity under inert annealing conditions. [𝐼] and [𝑉] are the local interstitial and vacancy concentrations, and [𝐼]∗ and [𝑉]∗ are the interstitial and vacancy concentrations under equilibrium conditions. 𝑓𝐼 and 𝑓𝑉 are the diffusion fractions mediated by interstitials and by vacancies respectively, and by definition 𝑓𝐼 + 𝑓𝑉 = 1. (2.5) These formulae are capable of catching the diffusion behaviors for both equilibrium and non-equilibrium diffusion. Point defect concentrations can be disturbed by using either external sources (such as ion implantation, oxidation, nitridation and proton irradiation) or internal sources (such as the incorporation of another species). 14 Before further discussions, we need to clarify some terminologies to avoid any confusions. Table 2.1 shows the definitions of several diffusion conditions. 𝑛𝑖 is the intrinsic carrier concentration, which is a material parameter and depends strongly on the temperature. In our study, the P concentrations are normally higher than 𝑛𝑖 , so the P diffusion of this work is extrinsic diffusion. Table 2.1 Terminology used for diffusion conditions. Terminology Conditions Intrinsic diffusion Doping level < 𝑛𝑖 Extrinsic diffusion Doping level > 𝑛𝑖 Unperturbed thermal Equilibrium [𝐼] = [𝐼]∗, [𝑉] = [𝑉]∗, [𝐼][𝑉] = [𝐼]∗[𝑉]∗ I and V recombine to a steady state [𝐼] ≠ [𝐼]∗, [𝑉] ≠ [𝑉]∗, [𝐼][𝑉] = [𝐼]∗[𝑉]∗ Non-steady state [𝐼][𝑉] ≠ [𝐼]∗[𝑉]∗ Historically, there have been numerous studies on the values of 𝑓𝐼 and 𝑓𝑉 for common dopants in Si. Fahey et al. provided comprehensive analysis to derive the boundary values of 𝑓𝐼 or 𝑓𝑉 for P, Sb and As [35]. Ural et al. offered improvements by carrying out numerical calculations and using SIMS rather than spreading resistance analysis which was used in earlier work [36]. Although there are debates on the values of 𝑓𝐼 and 𝑓𝑉 in some cases, such as in Si self-diffusion [27–29], this theoretical framework has been generally accepted. More details can be found in Ref.[23]. In the following sections, we will use the theoretical framework as discussed above. The determination of 𝑓𝐼 and 𝑓𝑉 for common dopants in Si has always been a central issue, as summarized by Pichler [23]. One important fact that we should always keep in mind is that, although these parameters should have specific values under certain conditions (e.g. at a certain temperature and for a certain dopant), they cannot be measured directly from experiments. All the 15 procedures that were used to indirectly obtain the 𝑓𝐼 and 𝑓𝑉 values involve experiments under perturbed conditions and subsequent analytical derivations. Therefore, the results obtained from various sources over the past several decades highly depend on the experiment conditions and the assumptions/approximations used. Without the capability to measure point defect concentrations directly, the best criterion that we can utilize so far is how well these 𝑓𝐼 and 𝑓𝑉 values agree with experimental data not only for P, but also for other dopants. 2.2.2 Fermi level effect in dopant diffusion Compared to diffusion in metals, there is a special effect in dopant diffusion in semiconductors, namely the Fermi effect, which comes from that fact that dopants exist as ionized impurities in semiconductors at room temperature and higher processing temperatures. For example, P exists in Si as 𝑃+ ions, and contributes one free electron per P atom. The Fermi effect is related to the diffusion mechanisms, where diffusion is carried out by both neutral and charged defects [22]. For example, interstitials in semiconductors can be neutral interstitials (𝐼0 ), or single negatively charged interstitials (𝐼− ), or double negatively charged interstitials (𝐼= ), or single positively charged interstitials (𝐼+ ), or double positively charged interstitials (𝐼++) etc. Same is true for vacancies. The concentration of each charged defect is closely related to the doping level. For example, when a semiconductor is n-type doped, there are extra free electrons, which promote the formation of negatively charged point defects. Then the diffusion of positively charged ions (𝑃+, for example) assisted by these negatively charged point defects will be enhanced consequently. The n-type dopants also reduce positively charged defect concentration. For p-type dopants, the effect is similar while all the charge states are just reversed. This effect is important when the doping level is larger than the intrinsic carrier concentration 𝑛𝑖 at the temperature of diffusion, and it is necessary to include the Fermi effect in this work. 16 The Fermi effect can be very strong, and introduce hundreds and thousands fold of increases in the dopant diffusivity for heavily doped semiconductors. Experimentally, if we examine the actual concentration profile of many diffused dopants, we find that the diffusion appears to be faster in the higher-concentration regions. A commonly accepted mathematical form to describe the Fermi effect on the diffusivity is polynomial. For example, for an n-type doped semiconductor, the diffusivity is expressed as: 𝐷𝑒𝑓𝑓 = 𝐷0 + 𝐷−(𝑛𝑛𝑖) + 𝐷=(𝑛𝑛𝑖)2, (2-6) where 𝑛 is the local free electron concentration. For P-type dopant, the individual diffusivity terms that are significant are 𝐷0, 𝐷+ and 𝐷++ and so forth. Each of these individual diffusivities can be written in an Arrhenius form as 𝐷 = 𝐷0exp (−𝐸𝐴𝑘𝐵𝑇) (2.6) with a preexponential factor 𝐷0 and an activation energy 𝐸𝐴 . Since this effect is thought to be related to the Fermi energy level that is set mainly by the doping level, it is often called the Fermi level effect. This effect has been well modeled and is included in TSUPREM-4TM [40]. In our modeling of coupled diffusion and segregation, we establish our own Matlab code, which includes this effect. 2.3 Methods of dopant profile engineering Dopant profile engineering during semiconductor fabrication can be categorized into two steps: one is the initial introduction of dopants such as pre-deposition, ion-implantation, in-situ doping etc., the other is the control of diffusion during high temperature process steps. In the case of PNP SiGe HBT devices where the phosphorus profile in base region is of central importance, the mission of the first step is to introduce the dopant with a high and narrow profile. For example, the atomic layer deposition (ALD) developed recently is such a technique for this purpose [41]. In this work, 17 the P introduction method is the in-situ doping during the epitaxial growth, which is the method used in industry. Our goal is to find out how to control P diffusion after the initial introduction step. A method to suppress dopant diffusion should satisfy two requirements simultaneously: it can reduce dopant diffusivity, and it should not bring significant negative changes to other material or device characteristics. Thanks to the nature that dopant diffusion “has to” rely on the assistance of either vacancies of self-interstitial atoms, it gives us opportunity to control diffusion. For an analogy, cars rely on bridges to cross a river, so by controlling the number of bridges it is possible to control the transport capacity. Since dopant diffusion relies on the help of point defects, we can control dopant diffusion by disturbing point defect concentrations. Point defect concentrations can be disturbed externally or internally as mentioned in the previous section. Excess vacancies or interstitials can be injected externally by gas ambient annealing (e.g. ammonia or oxygen [35]), or proton radiation [42]. However, proton radiation increases both interstitial and vacancy concentration, thus contradicted our technological mission to reduce P diffusion. Thermal nitridation is a possible method to suppress P diffusion in SiGe since it reduces interstitial concentrations, and P diffusion in Si and low Ge content SiGe alloys is mostly via interstitial mechanism. The incorporation of alien neutral impurities is also possible to suppress dopant diffusion. For example, carbon can suppress P and B diffusion in Si [43]. The current understanding is that carbon diffuses via interstitial mechanism (or kick-out mechanism), thus reducing the concentration of self-interstitials [11]. We will discuss this method and our experiments in more detail in Chapter 4. One thing to be noted is that, to our best knowledge, there has been no experimental study on the P diffusion mechanism in SiGe yet. The literature shows that P diffusion in Si is exclusively via interstitial while its diffusion in Ge is via vacancies. A natural logic may deduce that the P diffusion in SiGe should evolve from an interstitial mechanism to a vacancy mechanism with the increasing Ge molar fraction. In Chapter 4 and Chapter 5, we see that this hypothesis is supported by our quantitative experimental evidence. 18 Experiment and modeling methodologies 3.1 General procedure of diffusion study in this work To investigate diffusion, we designed, performed and compared diffusion behaviors of phosphorus in SiGe under various conditions. Some iterations were needed for the structure and diffusion thermal budget design. The standard procedure in this work is shown in Figure 3.1. Figure 3.1 The standard procedure of this work. 19 3.2 Sample structure design and growth 3.2.1 Strain and critical thickness for SiGe epitaxial heterostructures In a SiGe HBT, strained SiGe is used in the base region to create heterojunctions between Si on both sides. To preserve the crystallinity, strained SiGe layers are grown on top of Si in such a way that the lateral (in-plane) lattice spacing is maintained while the vertical (out-of-plane) lattice spacing can expand or contract freely. This can be done by epitaxial growth using a chemical vapor deposition (CVD) system, which is widely used in the semiconductor industry. Figure 3.2 The critical thickness vs. Ge fraction for pseudomorphic Si1−xGex layers grown on bulk (100) silicon . Figure reprinted from Ref.[1] with permission. Ge has a 4.2% larger lattice constant than Si, and the lattice constants for SiGe alloys obey Vegard’s law. The lattice mismatch between Si and SiGe impose a thickness limit, called the “critical thickness”, below which the strain can be maintained (see Figure 3.2). Beyond the critical 20 thickness, defects and misfit dislocations form to relieve the strain, which degrades the device performance. A higher difference in Ge content results in a larger mismatch strain, thus a smaller critical thickness. This practical thickness limit depends on the epitaxial growth condition, the thickness of Si capping layer on top and the subsequent thermal budget [44]. This metastable regime is the regime where most industry CMOS and HBT devices belong to. In this work, the SiGe layers were designed to be in the metastable regime. Therefore, we needed to design the structures and the thermal annealing conditions such that the strain in the SiGe layers is fully maintained. We also needed to monitor the strain status closely in the diffusion experiment. The design of SiGe layer thickness was based on two competing factors: it should be thick enough to accommodate sufficient P diffusion that SIMS can detect, and the layers can't be too thick to cause strain relaxation. A thicker layer can accommodate more diffusion, so it helps to obtain better accuracy in diffusivity extraction. Carbon incorporation is helpful to increase the practical thickness limit as it helps to reduce the compressive strain due to the much smaller carbon size in a SiGe lattice. The equilibrium critical thickness and the practical thickness limit both decrease rapidly with the increasing Ge fraction. Ge18C000 tends to relax the most as it has the highest Ge fraction in this study and no carbon. The Matthews-Blakeslee critical thickness model gives the equilibrium critical thickness, which is a very conservative thickness limit, below which the layers will not relax even with very high thermal budgets [45]. For 18% Ge, the Matthews-Blakeslee equilibrium critical thickness is less than 20 nm, which is too narrow to give sufficient room for dopant diffusion study in strained-SiGe. In practice, strain relaxation happens above a much thicker limit than the Matthews-Blakeslee critical thickness limit, as layers can be under a metastable state before strain relaxation starts. Therefore, some trials were necessary to find the thickness and thermal budget range without strain relaxation. Test structures with different SiGe layer thicknesses were grown and annealed. The strain statuses before and after anneals were checked by X-ray diffraction (XRD). Based on these trial 21 experiments, we found that 80nm 𝑆𝑖0.82𝐺𝑒0.18 can stay fully strained after diffusion experiments using the thermal budgets of our interest. There are two sets of samples studied in this work. One set of samples have uniform SiGe or SiGe:C regions where phosphorus diffusion happened. The other set of samples were used to study the segregation phenomenon across graded SiGe regions, thus the material composition in these region were designed to be non-uniform. For the latter set of samples, the maximum Ge molar fraction is 0.18. The thickness and thermal budget requirements for the latter set of samples are less stringent, as the averaged Ge concentration is only half of the peak Ge concentration. 3.2.2 Sample structure and growth All the designs are shown in more detail in the corresponding chapters. All the epitaxial layers were grown on 8 inch (100) Czochralski p-type Si wafers in an ASM Epsilon 2000 reduced-pressure chemical vapor deposition (RPCVD) reactor at Texas Instrument Deutschland GmbH in Germany. The epitaxial growth temperature of all structures was 550 ℃. 3.3 Annealing of diffusion samples The selection of the annealing tool was based on the thermal budget design (temperature and time). For technological relevance, the industry mainstream annealing tools are rapid thermal annealing (RTA) tools, which have fast ramp rate of about 100 ⁰C/sec. Soak RTA anneals normally last from one second to a few minutes at the peak temperature. However, if the annealing time at the peak temperature is too short, there will not be sufficient time for diffusivity extraction at this specific temperature. Therefore, all our samples were annealed for no less than 15 seconds. In terms of temperature choice, mainstream industry uses 1000 to 1100 ⁰C in the highest temperature anneals, so we used 1000 to 1050 ⁰C in the study of P diffusion in Chapter 4 and 5. In this work, three annealing tools were used. 22 1) In Chapter 4, inert anneals for carbon and Ge impact studies were performed in nitrogen ambient using an Applied Material Radiance rapid thermal processing tool at 1000 ⁰C for 15 to 120 seconds at Texas Instrument Deutschland GmbH in Germany . 2) In Chapter 5, ammonia anneals were performed using an AccuThermo AW610 atmospheric rapid thermal processing (RTP) system at Stanford Nanofabrication Facility (SNF), which has the rare capability of using NH3 as an ambient and can anneal wafer pieces. 3) In Chapter 6, inert anneals were performed using a Linkam Scientific Instruments’ high temperature stage system TS1200 in nitrogen ambient at 900 ⁰C in our lab at UBC. The reason for using the heating stage was to compare our results with literature data in similar temperature ranges and longer anneal times ( ≥10 mins), which the RTP tool at TI was not able to achieve. 3.4 Dopant profile characterization: secondary ion mass spectrometry The central task of this thesis is to investigate the P diffusion in strained SiGe. The thin thickness of strained SiGe layer, tens of nanometers, poses a requirement for precise dopant profiling on the scale of angstroms. Secondary ion mass spectrometry (SIMS) is the most appropriate analytical technique for this purpose, because it has the highest detection sensitivity (1017 cm-3) for measuring elemental concentration, and is able to profile in the depth dimension with sub-nanometer precision. SIMS is based on the fact that charged atomic and molecular species are ejected from the surface of a sample under heavy particle bombardment. Those secondary ions are then differentiated by a mass analyzer, which produce mass spectra indicating the elemental composition in the sample material. For diffusion study, the depth profiling of an interested element is crucial. Depth profiles in SIMS analysis can be obtained by monitoring the secondary ion count rate of selected elements as a function of time. To convert the time axis into depth, a SIMS analyst uses a profilometer to measure the sputter crater depth simultaneously. In this way, the elemental 23 composition as a function of depth is obtained, which constitute the profiles of elements. In this study, it is the phosphorus, carbon, silicon and germanium profiles which are obtained simultaneously for each sample. All the SIMS measurements were performed by Evans Analytical Group, which is the industry leader for commercial SIMS analysis. The samples were sputtered with 1 KeV Cs primary ion beam obliquely incident on the samples at 60° off the sample surface normal. The sputter rate was calibrated by a stylus profilometer measurement of the sputtered carter depth and corrected on a point-by-point basis. The measurement uncertainty for Ge fraction is ±1 at.%. For P and C, the concentration uncertainty is ±3% and ±10%, respectively. The depth/thickness uncertainty is approximately 5%. Generally, the P dose variation for each structure is less than 10% among as-grown and annealed samples. It should be noted that carbon concentrations measured from SIMS are total carbon concentrations regardless of the carbon positions in the lattice. 3.5 Strain characterization by X-ray diffraction As discussed earlier, what we studied was P diffusion in strained SiGe or SiGe:C. Long thermal anneals or thick strained layers may cause strain to relax. As strain impacts diffusion, it is important to monitor the strain status before and after diffusion closely. The most suitable characterization technique for strain measurements in our work is X-ray diffraction. 3.5.1 Principles of X-ray diffraction X-ray diffraction has been used to investigate semiconductor materials for decades. It has many unique advantages including its non-destructive nature, a good match of X-ray wavelength to the atomic scale of crystal structure and the rapid collection of statistically significant data [46]. A simple formulation of the “Bragg condition” is 2𝑑 sin𝜃 = 𝑛𝜆 (3.1) 24 where d is the spacing between diffracting planes, 𝜃 is the incident angle, n is any integer, and λ is the wavelength of the incident X-ray beam. This condition can be visualized as in Figure 3.4. Figure 3.3 A cross-sectional illustration of Bragg condition. When X-rays encounter atoms, a portion of its intensity is scattered by electrons, the scatterers. The scattered X-ray can be collected to generate a unique image called a “diffraction pattern”, which can be interpreted to reveal the atomic structure of the sample. For example, an intensity peak corresponds to a specific crystal plane inside the material and its orientation is featured by the Bragg angle. Therefore, the measurement of the angle of the X-ray beam is essential to the accuracy of this method. Modern high-resolution X-ray diffraction (HRXRD) can give results with the angular precision at tens of arc-seconds, which is used in this work. One important feature about XRD measurement is that it only gives us information about the crystal structure, not the identity of atoms. This is exactly the opposite of what we can obtain from SIMS analysis. 25 3.5.2 Reciprocal space mapping Reciprocal space mapping (RSM) is a special technique in XRD analysis. It generates a 2-dimensional cut of the 3-dimensional reciprocal space of the sample. In this way, the lattice strain status can be directly revealed in the result, which is usually an X-ray intensity contour. The principle of RSM can be illustrated as in Figure 3.4. Figure 3.4 Two dimensional reciprocal lattice and the Ewald sphere representation. In a reciprocal space map, the alignment of the diffraction spot from a SiGe layer with respect to the spot from the underlying Si substrate shows whether this layer is fully strained (i.e., with the same in-plane lattice constant) with the Si substrate. In a fully strained epitaxial layer, its in-plane lattice constant equals to that of the substrate. Thus, their reciprocal lattice points will align vertically in the reciprocal space map. If they are not aligned, it means that their in-plane lattice constants are different and relaxation has occurred. 26 In this work, XRD measurements were performed on the as-grown and the annealed samples to confirm the strain status and substitutional carbon molar fraction. All the measurements were performed using a PANalytical X’Pert PRO MRD with a triple axis configuration in Advanced Materials and Process Engineering Laboratory (AMPEL) at UBC. (113) reciprocal space mapping was performed to further confirm the strain status for the sample with the highest thermal budget. 𝜔 − 2𝜃 rocking curves of (004) Bragg reflection were measured to determine the substitutional carbon molar fraction for each sample. The X-ray tube was operated at 45 kV and 40 mA in the line focus mode. The Cu − Kα1 wavelength (λ = 1.5406 Å) was selected with a monochromator. 3.6 Computer modeling and simulations One of the greatest applications of a computer in scientific or engineering research is that it can simulate natural processes and make predictions. Its application in diffusion studies is an example. Simple analytical solutions such as complementary-error function or Gaussian function do not apply to more complicated diffusion problems and many industry problems like the case of HBTs. In this study, we needed to include the concentration-dependent diffusivity (i.e. the Fermi-effect), the impact of other impurities (carbon and germanium in this work), and dopant segregation. Therefore, the dopant diffusivity has to be calculated based on the local material composition and impurity profiles. In such cases, numerical analysis to simulate a natural process becomes a must. What it requires are mathematical models for these physical processes and the corresponding parameters, e.g. the Fick’s law and diffusivity in the case of diffusion. As long as the models and parameters are accurate, computer simulations can perform “virtual experiments”, which have many benefits to the industry R&D and design such as the reduction in cost, risk and time. On the other hand, computer simulations and fitting can also be part of the experiment analysis procedure, e.g. the extraction of diffusivity in our case. 27 3.6.1 Simulations with process simulator TSUPREM-4TM TSUPREM-4TM [40] is an industry mainstream two-dimensional (2-D) technology computer aided design (TCAD) tool to simulate processes used in the manufacturing of silicon integrated circuits. It simulates the incorporation and redistribution of impurities in a 2-D cross-section of a silicon wafer. As the diffusion and segregation in this work are one-dimensional in the direction normal to the wafer surface, a 2-D finite element analysis (FEA) tool is sufficient. A typical simulation structure consists of a number of material regions. Each material can be doped with multiple impurities. It also simulates the distribution of point defects in semiconductor regions during diffusion and their effects on the diffusion of impurities. The major reason that we used TSUPREM-4TM in this work is that it has well-calibrated models and material parameters based on several decades of industrial practices, such as the defect and diffusion models and the diffusivities of common dopants in Si and SiGe. In Chapter 4 and Chapter 5, the P diffusivity in each annealing condition was extracted by TSUPREM-4 fitting in terms of a ratio with regards to a reference diffusivity. TSUPREM-4 took the as-grown SIMS profiles as the pre-annealed profiles. The diffusivity used was the TSUPREM-4 default diffusivity of P in Si multiplied by a fitting parameter, i.e. the diffusivity ratio (𝑅𝐶𝐴,𝑆𝑖1−𝑥𝐺𝑒𝑥 ≡𝐷𝐴𝑆𝑖1−𝑥𝐺𝑒𝑥:𝐶𝑦𝐷𝐴𝑆𝑖1−𝑥𝐺𝑒𝑥 in Ch.4 and 𝑅𝑁𝑇𝐷 ≡𝐷𝑃𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝐷𝑃𝑖𝑛𝑒𝑟𝑡 in Ch.5). The fitting parameter value, i.e. 𝑅𝐶𝐴,𝑆𝑖1−𝑥𝐺𝑒𝑥 or 𝑅𝑁𝑇𝐷, was determined until the best match of simulated diffusion profile to the post-annealed SIMS profile was achieved. This way, we could study the impact of carbon or nitridation using the diffusivity ratio 𝑅𝐶𝐴,𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑅𝑁𝑇𝐷, respectively. TSUPREM-4TM has the capability to accommodate some user-defined models. However, it was not able to take the Ge concentration gradient introduced segregation flux into account. Therefore, we had to find another method to simulate the coupled diffusion-segregation 28 phenomenon. We used MatlabTM to numerically solve the diffusion and segregation transport equations using the finite difference method. 3.6.2 Finite difference method Equations describing mass transport phenomena are often partial differential equations, and we usually use the finite difference method for numerical analysis. The domain of interest is partitioned in space and in time and approximations of the solution are computed at the space or time points. This way, the continuous differential equations are solved by approximating them with “difference equations”, in which the derivatives are approximated by finite differences. Figure 3.5 Illustration of finite difference method for solving mass transport equations. Let’s consider a small spatial element (such as a box shape shown in Figure 3.5) where there are mass fluxes coming in and going out. The principle of mass conservation can be applied, which means “what goes in and doesn’t go out stays there”. For one-dimensional case, it can be written as ∆𝐶 ∙ ∆𝑥 = (𝐹𝑖𝑛 − 𝐹𝑜𝑢𝑡) ∙ ∆𝑡 , (3.1) where ∆𝐶 is the change of concentration due to the difference of in-flux 𝐹𝑖𝑛 and out-flux 𝐹𝑜𝑢𝑡 during time interval ∆𝑡. Using Fick’s first law as an example, the mass fluxes can be expressed as 𝐹𝑖𝑛 = −𝐷𝐶𝑖−𝐶𝑖−1Δ𝑥 and 𝐹𝑜𝑢𝑡 = −𝐷𝐶𝑖+1−𝐶𝑖Δ𝑥 (3.2) 29 where 𝐶𝑖−1 and 𝐶𝑖+1 are the concentrations of the adjacent spatial elements (see Figure 3.6). Thus, the Fick’s first law has also been expressed as difference equations. Figure 3.6 Illustration showing how diffusion flux such as Fick’s first law is expressed by finite difference equations. Combing Eq.(3.1) and Eq.(3.2), we can obtain a numerical expression of Fick’s law as 𝐶𝑖𝑡+∆𝑡 = 𝐶𝑖𝑡 +𝐷∆𝑡∆𝑥2(𝐶𝑖−1 − 2𝐶𝑖 + 𝐶𝑖+1) , (3.3) where the concentration of element 𝑖 at time 𝑡 + ∆𝑡 is explicitly expressed by the information at time 𝑡. By solving Eq.(3.3) across the desired space and time range, we can obtain the numerically solved diffusion profile. This procedure is illustrated in Figure 3.7. The method described here is sometimes called the “FTCS explicit method”: Forward in Time, Central in Space. The required condition for this method to be stable and can get convergent result is: ∆𝑡 ≤12∆𝑥2𝐷 (3.4) There are several other methods of finite difference analysis techniques, which will not be discussed here. 30 Figure 3.7 Illustrations of numerical calculation procedures across spatial and time domain. In this work, we used the finite difference method and implemented it in MATLAB to simulate the diffusion-segregation phenomenon in inhomogeneous materials. 31 The effect of carbon on phosphorus diffusion in SiGe As mentioned in 1.1.3, carbon incorporation has been widely used in NPN SiGe HBTs to control boron diffusion in SiGe, but there are fewer studies available on P diffusion in SiGe:C. In Chapter 4, we present a systematic and quantitative study of carbon impact on P diffusion in strained SiGe through experiments and modeling. This topic is of great industry relevance, which can give guidance to the design of SiGe:C composition, structure and thermal processing. 4.1 Introduction The fact that carbon can retard B and P diffusion in Si can be explained by two points: 1) carbon can reduce interstitial concentrations of Si and 2) B and P have very strong preferences on the interstitial mechanism, which can be represented by the value of 𝑓𝐼, i. e. 𝑓𝐼𝐵,𝑆𝑖 ≈ 1, 𝑓𝐼𝑃,𝑆𝑖 ≈ 1. For other common dopants, we have 𝑓𝐼𝐴𝑠,𝑆𝑖 ≈ 0.4 and 𝑓𝐼𝑆𝑏,𝑆𝑖 ≈ 0.02 [22]. The value of 𝑓𝐼 is crucially important from a technological point of view because it indicates whether and how we can change the diffusivity of a specific dopant. The incorporation of an extra chemical element is a typical technological choice for the control, mostly the suppression, of dopant diffusion if that element can change the point defect concentrations in a desired way. Besides the impact on diffusion, two extra properties are desired: it is electrically inactive (neutral impurity) and its impact on material properties is negligible at the given concentration. Based on these criterions, carbon became the ideal solution to suppress interstitial diffusers, such as B and P, as carbon can reduce the concentration of interstitials. Also, carbon is a neutral element for Si/SiGe system so that it does not affect charge distribution. The size of carbon atom is smaller than Si and Ge so that it slightly compensates the Ge-induced lattice-mismatch strain and increases the stability and critical thickness of strained SiGe layer [47]. 32 It has been reported that carbon of 0.1 at.% concentration can suppress B and P diffusion in Si by a factor of 10 and even more with higher concentration [48]. When this method is extended from Si to SiGe, the situation has changed. The carbon incorporation method is still effective for boron, but not for phosphorus. This drew us to question the underlying reason for their different responses to the same stimulus. To address this question, let’s begin with literature review on the effect of carbon on dopant diffusion. 4.2 Literature review: the effect of carbon on dopant diffusion in Si and SiGe 4.2.1 Carbon in Si and SiGe Carbon is a very common impurity in silicon wafers. It can be introduced either from the polysilicon used in crystal growth, or through contamination during this process [23]. Typical carbon concentrations in silicon wafers used in VLSI technology are now below 1016𝑐𝑚−3. When carbon is used for bandgap, stress or diffusion engineering, epitaxial deposition or solid-phase epitaxial regrowth of carbon-doped layers allow the incorporation of carbon up to 1021𝑐𝑚−3 (2 at.%) in Si and SiGe [49], [50]. Carbon has many formats of existence in a silicon lattice. It can exist as substitutional carbon (𝐶𝑆) or interstitial carbon (𝐶𝐼) in the silicon lattice, or, it can form small or large complexes with other atoms or carbon clusters. Newman and Wakefield concluded that, in bulk monocrystalline silicon, carbon resides predominantly on substitutional sites based on the observations that the diffusion coefficient and the activation energy of carbon are similar to those of the dopants [51]. The substitutional character of carbon in silicon was experimentally demonstrated by Baker et al. [52] and Windisch and Becker [53] who found that the presence of carbon causes a reduction of the lattice parameter to an extent close to what can be expected from a linear interpolation between pure silicon and 𝛽-SiC. A second 33 prominent atomic configuration is that carbon atoms in the interstitial sites [54], which usually results from the reaction of substitutional carbon with a silicon self-interstitial. Because of the difference in covalent bond radius between Si and C atoms, C atoms in a Si lattice are surrounded by a strain field, with substitutional carbon atoms under tensile strain and interstitial ones under compressive strain fields [55]. The minimum formation energy of an interstitial carbon in split (100) configuration is 4.6 eV, which exceeds the formation energy of the substitutional carbon by 3 eV [56], which explains why carbon atoms are favored to occupy substitutional sites. In terms of epitaxially grown Si:C and SiGe:C films, the ratio of substitutional carbon concentration with respect to total carbon concentration (𝐶𝑠𝐶𝑡𝑜𝑡𝑎𝑙) depends on the growth methods and conditions, and the value of this ratio was experimentally measured in this work. 4.2.2 Substitutional carbon incorporation in epitaxial grown Si1-yCy and Si1-x-yGexCy Compared to Si and SiGe epitaxy without carbon, epitaxial growth of C-containing Si1-yCy or Si1-x-yGexCy alloys requires additional problems to be addressed. Unlike the Si-Ge system, silicon and carbon are only miscible to a very small degree under equilibrium conditions [6]. According to the Si-C phase diagram, stoichiometric SiC (silicon carbide) is the only stable compound [57]. Any alloy with a smaller C concentration is thermodynamically metastable. Such alloy layers can be obtained by kinetically dominated growth methods, such as molecular beam epitaxy (MBE), [58], solid-phase epitaxy [59], or CVD [60]. All of these methods generally work under far from thermodynamic equilibrium conditions, allowing a kinetic stabilization of metastable phases. Although the bulk solubility of carbon in silicon is small, 3 × 1017𝑐𝑚−3, at the melting point of silicon [61], epitaxial layers with more than 1%, 5 × 1020𝑐𝑚−3,carbon can be achieved. However, during an epitaxial growth step, it is not the equilibrium bulk solubility that is important but the “surface solubility” [62]. 34 Substitutional carbon fraction or substitutionality ( 𝐶𝑠𝐶𝑡𝑜𝑡𝑎𝑙) is strongly influenced by the growth conditions [50], [63]. Similar trends were observed for MBE and CVD growth, that lower growth temperature with a higher growth rate resulted in higher substitutional carbon concentration. Figure 4.1 shows an example of the relation between carbon substitutionality and growth conditions with a low-pressure rapid thermal chemical vapor deposition tool [60]. The carbon substitutionality was determined by the cross-comparison between the XRD and SIMS measured carbon concentration. The XRD results represent substitutional carbon concentration while SIMS results represent total carbon concentration. The diagonal line represents the condition of fully substitutional incorporation (i.e. 𝐶𝑠𝐶𝑡𝑜𝑡𝑎𝑙= 100%). It is clearly shown that decreasing the growth temperature to 550 ℃ dramatically improves the carbon substitutional incorporation. In this work, due to the presence of Ge, and the 1% SIMS error in Ge molar fraction measurement, we can not use the cross-comparison between the XRD and SIMS measured carbon to determine the carbon substitutionality. Instead, we used the comparison between XRD of SiGe and SiGe:C and SIMS results to measure the carbon substitutionality, as discussed in 4.3.2. 35 Figure 4.1 Summary of the carbon substitutionality in epitaxial Si1-yCy for different growth temperatures and silane partial pressures. Figure reprinted from Ref.[60] with permission. Figure 4.2 shows the carbon substitutionality in epitaxial Si1-x-yGexCy for different growth temperatures and carbon sources. This trend is similar to the case of Si1-yCy. In this work, all our wafers were grown at 550 ℃ and the carbon concentration is up to 0.32%, which is well inside the region of complete substitutional incorporation. These carbon concentrations are comparable to the dopant concentrations. To distinguish these materials from the Si1-yCy or Si1-x-yGexCy systems discussed above, we will use the same notation as for dopants, that is Si:C or SiGe:C. This notation also helps to convey our assumption that the carbon concentration is too low to significantly affect band alignment and strain of the Si:C or SiGe:C layers [62]. 36 Figure 4.2 Summary of carbon substitutionality in epitaxial Si1-x-yGexCy for different growth temperatures and carbon sources. The Si1-x-yGexCy samples include Ge contents from 10 to 36%. Figure reprinted from Ref.[63] with permission. 4.2.3 Mechanisms of carbon in defect and diffusion engineering The physical mechanisms of the carbon impact on dopant diffusion lie in the undersaturation of Si self-interstitials [11]–[13]. The diffusion of carbon is generally assumed to proceed via an interstitial-assisted mechanism. Evidence for such a mechanism resulted from infrared-absorption experiments where it was observed that the production of vacancy-oxygen pairs and divacancies is increased in the presence of carbon [64]. A significant fraction of carbon diffusion via interstitial mechanism is also evidenced by the experiments of Kalejs et al. who observed enhanced carbon diffusion under conditions where the concentration of self-interstitials is enhanced by thermal oxidation or diffusion of phosphorus [65], [66]. Therefore, when a high and non-uniform 37 concentration of 𝐶𝑠 is grown in, the immobile 𝐶𝑠 needs to change to mobile 𝐶𝐼 for diffusion to happen [67]. The reaction of this step and its reverse reaction can be described as: 𝐶𝑆 + 𝐼 ↔ 𝐶𝐼 (4.1) 𝐶𝐼 + 𝑉 ↔ 𝐶𝑆 (4.2) The first reaction describes that a 𝐶𝑠 kicks a Si interstitial out from its site, and moves to the interstitial lattice site, also known as the kick-out mechanism. It consumes one Si interstitial per 𝐶𝑠 atom. It is worth noting that the transition from 𝐶𝑠 to 𝐶𝐼, and thus the consumption of Si interstitials, can happen without carbon diffusion. Rücker et al. studied the impact of C on B diffusion where no supersaturation of interstitials due to ion-implantation or other processing steps was present [12], [13]. The diffusivity of B was found to be reduced by about a factor of twenty for temperatures between 750 ℃ and 900 ℃ due to the presence of 1020𝑐𝑚−3 substitutional C. This effect has been attributed to a reduced density of Si self-interstitials [11], [68], [69]. These results show that carbon can be used as a diffusion suppressing agent if the target dopant is an interstitial diffuser. The use of C-rich layers to reduce B diffusion has been widely used in SiGe HBTs [70], [71]. Transistors with improved static and dynamic performance can be fabricated with epitaxial SiGe:C layers [72]. 4.2.4 Carbon impacts B and P diffusion in Si similarly It has been generally agreed that B and P diffusion in Si is almost exclusively via the interstitial-assisted mechanism, i.e. the interstitial-assisted diffusion fraction 𝑓𝐼 ≈ 1 [25, 38]. It means that any method which can reduce the interstitial concentration will suppress B and P diffusion in Si. According to these studies, the degree of the carbon suppression effect on B and P diffusion in Si is very similar, presumably because the value of 𝑓𝐼 for B and P are almost identical. 38 4.2.5 Carbon impacts B and P diffusion in SiGe differently A handful of studies are available on B and P diffusion in low Ge content SiGe. Zangenberg et al. and Christensen et al. both assumed that B and P were interstitial diffusers in SiGe [73], [74]. When the diffusion media changes from Si to strained-SiGe, the approach of using carbon to control diffusion has also worked quite well for B in SiGe bases for NPN HBTs [75]. It has been reported that carbon can retard B diffusion very similarly in strained SiGe with Ge molar fraction up to 20% [15]. However, the case for P is quite different. In 𝑆𝑖0.8𝐺𝑒0.2, the extent of suppression effect on P diffusion is similar for 0.06% (3 × 1019 cm-3) and 0.2% (1 × 1020 cm-3) carbon molar fractions, while for B diffusion, the higher 0.2% carbon molar fraction clearly retards B stronger than the 0.06% carbon [9]. The above studies suggest that the diffusion mechanisms of B and P in Si may be similar, but they are different enough in strained-SiGe to cause the difference of carbon retardation effectiveness. So far, there has been no systematic and quantitative study on the carbon retardation effect on P diffusion in SiGe, which is addressed in this work. Next, we present the experiment and modelling work on phosphorus diffusion in SiGe and SiGe:C. 4.3 Diffusion experiments and results 4.3.1 Sample structures and annealing conditions design To understand the P diffusion behavior and mechanisms in 𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦, a wafer matrix with various Ge and C contents was designed and grown by epitaxy. Figure 4.3 shows the epitaxial structure design, where 𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦 are sandwiched between two Si layers. To be relevant to PNP SiGe HBT applications, the x value was chosen to be 0.18 and the y values ranged from 0 to 0.0032. The naming convention of the structures is to use the Ge and C 39 percentages to name a 𝑆𝑖/𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦/𝑆𝑖 structure. For example, Ge18C018 is used to name a 𝑆𝑖/𝑆𝑖0.82𝐺𝑒0.18: 𝐶0.0018/𝑆𝑖 structure, which has 0.18% carbon as measured by SIMS. For each type of structure, several wafers were grown under identical growth conditions to be used as as-grown wafers and annealed wafers. Inert anneals were performed in nitrogen ambient using an Applied Material Radiance rapid thermal processing tool at 1000 ⁰C for 15 to 120 seconds. It should be noted that the concentration of substitutional C in epitaxial layers can exceed the solid solubility limit (e.g. 3.5 × 1017 cm−3 near the melting point) by several orders of magnitude [76]. Under optimized growth conditions, almost all C atoms can be incorporated on substitutional sites, which are the lowest energy configuration for isolated C atoms [50]. In this study, we confirmed this condition by X-ray diffraction (XRD) measurements for all as-grown and annealed samples. The SiGe layers were grown pseudomorphically on Si substrates so that biaxial compressive strain was maintained. Phosphorus was introduced in-situ by epitaxial growth with peak concentrations in the range of 2~4 × 1019 cm−3 . The choice of SiGe layer thickness was based on two considerations: it should be thick enough to accommodate sufficient P diffusion within SiGe region for better accuracy in diffusivity extraction while the strain should be maintained without relaxation. It is known that the equilibrium critical thickness reduces rapidly with increasing Ge fraction [44]. Ge18C000 is the one that tends to relax the most according to the critical thickness criteria. For 18% Ge, the well-known Matthews-Blakeslee model shows that equilibrium critical thickness is less than 20 nm, which is too narrow to give sufficient room for a dopant diffusion study in strained-SiGe. In practice, strain relaxation happens above a much thicker limit, which depends on epitaxial growth condition, the thickness of Si capping layer on top and the annealing thermal budget. To find out the thickness range without strain relaxation, test structures with different SiGe layer thicknesses were grown and annealed. The strain statuses before and after 40 anneals were checked by XRD. Based on these trial experiments, we found that 80nm 𝑆𝑖0.82𝐺𝑒0.18 can stay fully strained after diffusion experiments using the thermal budgets of our interest. Figure 4.3 Schematic diagram of the structure design and splits in C mole fractions. The dotted line indicates the original wafer surface. The carbon concentration split is from 0.06 to 0.32 at.%, which corresponds to 3e19 to 1.6e20 cm-3. The Ge and C profiles are designed to be constant with a certain thickness as shown above. 4.3.2 Strain and depth profile analysis (113) reciprocal space mapping (RSM) was performed on as-grown and annealed samples at room temperature. Figure 4.4 shows the RSMs of some annealed samples with the largest thermal budget for each type of structure. These results clearly show that the in-plane lattice constants of the 𝑆𝑖0.82𝐺𝑒0.18 and 𝑆𝑖0.82𝐺𝑒0.18: 𝐶𝑦 layer are the same as that of the Si substrate after annealing, which means that these two layers are fully strained to the substrate during diffusion. The (113) RSMs were done on all structures and annealing conditions, and no strain relaxation was observed 41 for the structures after diffusion in this work. These results confirm that P diffusion happened in fully strained SiGe or fully strained SiGe:C. In the following discussion, SiGe stands for fully strained SiGe unless otherwise noted. Figure 4.4 (113) reciprocal space mapping of annealed samples with the largest thermal budget for each type of structures. a) Ge18C000; b) Ge18C009; c) Ge18C018; d) Ge18C032. 𝜔 − 2𝜃 rocking scans of (004) Bragg reflection were performed to confirm the substitutional carbon fractions for each structure. Since XRD measurements can not measure interstitial carbon, 42 for a SiGe:C sample, (004) peak shifts only depend on the Ge molar fraction and the substitutional carbon molar fraction CS. In theory, the Ge molar fraction can be obtained using the SIMS data; However, the ±1 atom % uncertainty in SIMS Ge molar fraction can lead to a large CS uncertainty. Instead, we extracted the as-grown CS by comparing the XRD peak difference between an as-grown SiGe control sample and an as-grown SiGe:C sample, which were grown with the same epitaxy recipe except for carbon incorporation. To eliminate the influence of the substrate wafer miscut and the tilt of the epitaxial films, XRD scans were performed with the wafer oriented at 𝑝ℎ𝑖 = 0° and 180° . By matching the measurements with simulations using the PANalytical Epitaxy software package, the substitutional carbon molar fractions were extracted. By comparing the XRD simulations and SIMS data, we concluded that carbon in as-grown samples is almost exclusively substitutional. Comparing the as-grown and annealed samples, we found that the substitutional carbon molar fraction is only reduced by less than 0.005% for Ge18C020 samples (see Figure 4.5). That means the substitutional carbon in annealed samples is more than 97% of that in as-grown samples. In fact, our growth condition has been optimized to achieve complete substitutional incorporation [62], [63]. Therefore, we approximate that the substitutional carbon 𝐶𝑠 molar fraction is 100% in the modeling section, i.e. 𝐶𝑠 = 𝐶𝑦. 43 Figure 4.5 XRD (004) scan of SiGe and SiGe:C as-grown samples showing carbon molar fraction change. SIMS measurements were performed for all as-grown and annealed samples. Figure 4.6 shows an example of SIMS profiles for as-grown and annealed wafers. Both the as-grown and annealed P peaks are within the constant Ge and C region. The shape of the diffused P peak is close to a Gaussian distribution, indicating that the concentration dependence of P diffusivity in SiGe is not strong. There is about 1% Ge concentration difference between the as-grown and annealed samples, which is within the SIMS measurement error. Ge profile also broadens slightly due to Si-Ge interdiffusion. C diffusion is seen to be much faster than P and Ge diffusion. The Ge fraction measured by XRD is 17.9% ± 0.1%. The Ge fraction measured by SIMS is consistent with that obtained from XRD within measurement errors. 1.E+011.E+021.E+031.E+041.E+051.E+061.E+07-2400 -1400 -400 600 1600Intensity (cps)Relative angle (arcsecond)Ge18C000Ge18C006Ge18C009Ge18C018Ge18C032 44 Figure 4.6 As-grown and annealed Ge, P and C SIMS profiles of Ge18C018 samples. 4.3.3 P diffusivity ratio extraction The P diffusivity in each annealing condition was extracted by TSUPREM-4 fitting. TSUPREM-4 took the as-grown SIMS profiles as pre-annealed profiles. The diffusivity used was the TSUPREM-4 default diffusivity of P in Si, 𝐷𝑃𝑆𝑖, multiplied by a fitting parameter. The fitting parameter value was modified until a good match to the post-annealed SIMS profile was achieved. At the annealing temperature of 1000 ⁰C, our SIMS and fitting results showed that P diffusivity in 𝑆𝑖0.82𝐺𝑒0.18, 𝐷𝑃𝑆𝑖0.82𝐺𝑒0.18, is very close to the diffusivity of P in Si in TSUPREM-4, 𝐷𝑃𝑆𝑖, which is: 𝐷𝑃𝑆𝑖0.82𝐺𝑒0.18 ≈ 𝐷𝑃𝑆𝑖 (4.3) The default P diffusivity in TSUPREM-4 has been calibrated based on large amount of experimental data. We define 𝑅𝐶𝐴, 𝑆𝑖1−𝑥𝐺𝑒𝑥, the carbon impact factor, as the ratio of dopant element A diffusivity in 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦 over that in 𝑆𝑖1−𝑥𝐺𝑒𝑥: 101610171018101910200510152025300.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2P As-grownP AnnealedC As-grownC AnnealedGe As-grownGe AnnealedP, C Concentration (cm-3) Ge Molar Fraction (%)Depth (um)1000oC 70 sec 45 𝑅𝐶𝐴,𝑆𝑖1−𝑥𝐺𝑒𝑥 ≡𝐷𝐴𝑆𝑖1−𝑥𝐺𝑒𝑥:𝐶𝑦𝐷𝐴𝑆𝑖1−𝑥𝐺𝑒𝑥 (4.4) For each carbon molar fraction and each annealed sample, there is a specific value of 𝑅𝐶𝑃 extracted from SIMS profiles by TSUPREM-4 fitting. Figure 4.7 shows an example of the best fitting to a SIMS profile for 𝑅𝐶 extraction. The extracted 𝑅𝐶 values and the corresponding carbon molar fraction 𝐶𝑦 are plotted in Figure 4.8 in reference to the available literature 𝑅𝐶𝐵,𝑆𝑖1−𝑥𝐺𝑒𝑥 values(x ≤ 0.2). It can be seen that 𝑅𝐶𝑃,𝑆𝑖0.82𝐺𝑒0.18 are within the range of 0.43 to 0.85, which shows that carbon is still effective in retarding P diffusion in strained SiGe. However, the strongest retardation is in C molar fraction range of 0.05 to 0.1%, above which the retardation is less effective. This trend is different from P in Si, B in Si and B in 𝑆𝑖1−𝑥𝐺𝑒𝑥 ( x ≤ 0.2), where higher carbon gives more retardation. For C molar fraction at 0.20%, 𝑅𝐶𝑃,𝑆𝑖, 𝑅𝐶𝐵,𝑆𝑖 and 𝑅𝐶𝐵,𝑆𝑖1−𝑥𝐺𝑒𝑥 (x ≤ 0.2) are all around 0.1 or below, which shows much stronger retardation than C effect on P in 𝑆𝑖0.82𝐺𝑒0.18: 𝐶𝑦. This difference between 𝑅𝐶𝑃,𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑅𝐶𝐵,𝑆𝑖1−𝑥𝐺𝑒𝑥 clearly indicates that the diffusion mechanism for B and P in SiGe may be different. 46 Figure 4.7 Example of TSUPREM-4 fitting with SIMS profiles of as-grown and annealed samples Ge18C018. 00.511.50 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4P, 18% Ge, 1000 oC (This work)P & B, Si, 900 oC (Ref.)B, Si, 750-900 oC (Ref.)B, Si, 900 oC (Ref.)B, 10%Ge, 900 oC (Ref.)B, 20%Ge, 900 oC (Ref.)B, Si & 11%Ge, 1000 oC(Ref.)P, Si, 1000 oC (This work)RCCarbon Concentration (Atom %) Figure 4.8 Summary of the carbon impact factor for P and B diffusivity in Si1-xGex:Cy. Literature data are from Ref.[12], [15], [16], [48]. 101610171018101910200510152025300.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2P As-grownP AnnealedP FittingC As-grownC AnnealedGe As-grownGe AnnealedP, C Concentration (cm-3) Ge Molar Fraction (%)Depth (um)1000oC 70 sec 47 In the following, the difference between B and P diffusion will be discussed based on point-defect mediated diffusion mechanisms, i.e., interstitial-assisted and vacancy-assisted mechanisms. We will model the carbon effect, and extract the diffusion fractions mediated by interstitials and by vacancies, which are 𝑓𝐼 and 𝑓𝑉 respectively. With the carbon effect model and 𝑓𝐼 and 𝑓𝑉 values, the difference in P and B diffusion behaviors can be explained. Considering the relatively low concentration of carbon, we will treat Si and SiGe as the diffusion media material, as the role of carbon is to change the point-defect concentrations and stress in SiGe. 4.4 Modeling of carbon effect on P diffusion The effect of carbon on point defects is a dynamic process, and can be modelled by solving carbon-point defect reaction and transport equations [33, 45, 46, 55]. In major process simulators such as TSUPREM-4 [40] and Sentaurus Process [78], carbon effect is modeled using multiple reaction equations, such as Eq. (4.1) and (4.2). Higher order of carbon-defect complexes, typically involving more carbon and interstitial numbers, are sometimes used to model carbon effects [79]. Parameters such as chemical reaction rates, equilibrium concentrations and diffusivities of all mobile species are required to solve these equations. In principle, such an approach can describe the time evolution and space distributions of the point defects, and thus should be applicable to wide ranges of structures and annealing conditions. However, since these parameters are vastly beyond experimental verification, it is very difficult to calibrate these equations and use them in a reliable manner. Considering that in the current HBT industry, the SiGe:C structures and the annealing conditions are in a certain range, it is then meaningful to establish empirical models to capture the carbon effect in simpler formats, which are easier to calibrate and more time-efficient for 48 simulations. Therefore, our modeling will be based on the assumption that the local point defect concentrations recombine to a steady-state condition, which is [𝐼][𝑉] = [𝐼]∗[𝑉]∗ (4.5) Furthermore, it is also assumed that [𝐴𝐼]/[𝐴𝐼]∗ = [𝐼]/[𝐼]∗ and [𝐴𝑉]/[𝐴𝑉]∗ = [𝑉]/[𝑉]∗ (4.6) where [𝐴𝐼] and [𝐴𝑉] are the dopant-interstitial pair and dopant-vacancy pair concentrations, and * denotes their counterparts under equilibrium condition. Eq. (4.5) shows the point defect dynamics under perturbations: when the interstitial population of a material is decreased by X times, the vacancy population increases by X times. A similar approach has been applied to explain the effect of carbon on boron diffusion in Si and SiGe where B was treated as an exclusive interstitial-diffuser, i.e. 𝑓𝐼𝐵 = 1 [77]. Although it has been generally accepted that the carbon affects dopant diffusion through the interaction with point defects, a quantitative relation between carbon concentration and intrinsic point defect concentration is still lacking. Several studies on the modeling of carbon effects have been reported. Rucker et al. considered the "kick-out" reaction involving 𝐶𝑆, 𝐼 and 𝐶𝐼 (Eq.4.1) and expressed the equilibrium self-interstitial concentrations with regard to substitutional carbon concentration using a simple analytical function [75]: [𝐼]̃ [𝐼]∗≈11+𝑘𝑝 [𝐶𝑠] (4.7) where 𝑘𝑝 is the equilibrium constant for the reaction shown in Eq.(4.1). Sibaja-Hernandez et al. modeled the carbon effect considering carbon clustering and the corresponding reaction dynamics [77]. This approach requires more than 15 fitting parameters to solve the model system. In this work, the focus is to obtain an empirical and reliable estimation on the extent of interstitial suppression caused by carbon incorporation. Based on the trend of 𝑅𝐶𝐵,𝑆𝑖 and 𝑅𝐶𝐵,𝑆𝑖1−𝑥𝐺𝑒𝑥 data in Figure 4.8 and Eq. (4.5), we established the carbon impact on the ratios of the 49 time-averaged point defect concentrations ([𝐼]̃ and [𝑉]̃ ) over the corresponding equilibrium values ( [𝐼]∗ and [𝑉]∗) as: [𝐼]̃ [𝐼]∗=11+𝑘 [𝐶𝑠] (4.8) [𝑉]̃ [𝑉]∗= 1 + 𝑘 [𝐶𝑠] , (4.9) where k is a fitting parameter. In Eq. (4.8) and (4.9), the constant k is in front of the carbon molar fraction and it shows the strength of carbon in reducing interstitial concentration, and should not depend on dopant species. It effectively shrinks or expands the carbon molar fraction axis in Figure 4.8. In the following discussions, carbon is in molar fractions so that k is dimensionless. Figure 4.8 shows the carbon impact on B diffusivity in the temperature range of 750~1000 ⁰C and up to xGe = 0.2. It can be seen that the dependence of 𝑅𝐶𝐵 on the temperature and Ge fraction is negligible within the given ranges. Therefore, in this work, we assume that k is constant for xGe ≤ 20% within the range of 750~1000 ⁰C. Now we can combine Eq. (2.4), (2.5), (4.8) and (4.9) and describe the carbon effect on P diffusivity as 𝐷𝐴𝑒𝑓𝑓𝐷𝐴∗ = {𝑓𝐼𝐴 ∙11+𝑘 [𝐶𝑠] + (1 − 𝑓𝐼𝐴) ∙ (1 + 𝑘 [𝐶𝑠])} (4.10) Since the majority of the carbon is located in substitutional sites according to XRD results and the post-annealed carbon profiles remain as box-like with wide plateaus, we can replace the substitutional carbon fraction [𝐶𝑠] by total carbon fraction [𝐶𝑦] measured by SIMS. Therefore the Eq. (4.10) is translated to 𝐷𝐴𝑆𝑖1−𝑥𝐺𝑒𝑥:𝐶𝑦𝐷𝐴𝑆𝑖1−𝑥𝐺𝑒𝑥≡ 𝑅𝐶𝐴,𝑆𝑖1−𝑥𝐺𝑒𝑥 = 𝑓𝐼𝐴 ∙11+𝑘 [𝐶𝑦] + (1 − 𝑓𝐼𝐴) ∙ (1 + 𝑘 [𝐶𝑦]) (4.11) In this equation, both k and 𝑓𝐼 determine the shape of 𝑅𝐶(𝐶𝑦). In particular, the carbon dependence is strongly influenced by 𝑓𝐼, when 𝑓𝐼 is close to 1, see Figure 4.9. This model can 50 explain the “U-shape” of 𝑅𝐶 which we observed for P in SiGe. Figure 4.9 shows the relation between 𝑅𝐶 and self-interstitial suppression [𝐼]/[𝐼]∗ . When 𝑓𝐼 = 1 , carbon retards diffusion monotonically with the increase of carbon content. When 𝑓𝐼 = 0.99, the carbon retardation effect will saturate at [𝐼]/[𝐼]∗ ≈ 10 because the vacancy supersaturation will become dominant beyond this point. This feature implies that carbon can be used as an effective tool to “magnify” the small vacancy contribution. Figure 4.9 The relation between the interstitial suppression ratio ([I]/[I]*) and dopant diffusivity suppression ratio (RC = D/D*) for different diffusion mechanism under the assumption of [I][V] = [I]*[V]*. The above modeling only included the carbon impact on dopant diffusion in terms of changing the point defect concentrations. Moreover, carbon also affects the diffusion media strain by compensating compressive strain in SiGe and thus affects the equilibrium diffusivities of a dopant. The highest carbon molar fraction in this study is 0.32%, which means that the strain compensation 0.010.11101000.0100.1001.000RC[I]/[I]*𝒇𝑰 =00.40.60.900.950.980.991𝒇I = 0𝒇I = 1 51 introduced by carbon is equivalent to 3.2% of Ge, or about 0.13% in strain. However, these two types of carbon effects cannot be dissociated experimentally since they are always in effect simultaneously. Figure 4.10 shows a simplified illustration of the lattice configuration for typical structures grown on a Si substrate under different strain statuses. It can be noticed that the strained SiGe:C lattice is different from the relaxed Si:C or the relaxed SiGe or the strained SiGe lattices, which were previously studied. Therefore we can only estimate the strain compensation effect, commonly expressed as a change in the diffusion activation energy, based on previous studies on strain impact on dopant diffusion in strained-SiGe and relaxed-SiGe structures [80]. Figure 4.10 Simplified illustration of cross-section view of Si1-xGex and Si1-xGex:Cy lattice spacing under different strain status: (a) Si and fully strained Si1-xGex:Cy with x/y ≈ 10, and the out-of-place lattice constant is the same as the Si out-of-plane lattice constant, (b) fully strained Si1-xGex:Cy with x/y > 10 , (c) fully strained Si1-xGex and (d) fully relaxed Si1-xGex . Now we add an apparent strain factor to Eq. (4.11) and revise the model to: 52 𝐷(𝑃,𝑆𝑖𝐺𝑒:𝐶)𝐷(𝑃,𝑆𝑖𝐺𝑒)= {𝑓𝐼 ∙11+𝑘 [𝐶𝑦] + (1 − 𝑓𝐼) ∙ (1 + 𝑘 [𝐶𝑦])} ∙ 𝑒𝑥𝑝(−𝑞′∙𝜀𝑐𝑘𝐵𝑇) (4.12) where 𝜀𝑐 is the change of strain due to substitutional carbon ( 𝜀𝑐 < 0) , 𝑞′is the strain derivative with a unit of eV per unit strain, and 𝑘𝐵is the Boltzmann constant. Christensen et al reported an apparent activation energy of 𝑞′ = −13 𝑒𝑉 per unit strain, assuming that the strain 𝜀𝑐 is proportional to the Ge content by 𝜀𝑐 = −0.0418𝑥 [74]. Based on this value, 0.2% carbon will increase P diffusivity by about 10% due to the strain compensation effect. At 0.2% carbon, it was shown that the effective diffusivities of B and P were suppressed by more than a factor of 10, as shown in Figure 4.8. Therefore, in Eq. (4.12), the term involving point defects in the curly brackets is dominant. We will neglect the strain compensation effect of carbon (the exponential term in Eq. (4.12)) in the following sections, and use Eq. (4.11) in the rest of the work. After the model is established, we need to extract the parameters k and 𝑓𝐼 in Eq. (4.11). 4.5 Model parameter extraction Let’s consider the case of Si first. 𝑓𝐼𝑃,𝑆𝑖 ≈ 𝑓𝐼𝐵,𝑆𝑖 ≈ 1 has been generally agreed upon [17, 26]. Taking this assumption, we can extract 𝑘 from the 𝑅𝐶𝑃,𝑆𝑖 and 𝑅𝐶𝐵,𝑆𝑖 data, which overlay within measurement errors. Therefore, 𝑘 is independent of dopant species as discussed in section 4.2. Based on the fitting to this data (shown as the blue solid line in Figure 4.11), we estimate that 𝑘𝑆𝑖 ≈5500 for both B and P in Si. For SiGe, the literature data for boron shows that the carbon impact factor RC is very close to that in Si, i.e. 𝑅𝐶𝐵,𝑆𝑖1−𝑥𝐺𝑒𝑥 ≈ 𝑅𝐶𝐵,𝑆𝑖 for 0 ≤ 𝑥 ≤ 0.2 with the C range up to 0.2%, which indicates that 𝑓𝐼𝐵,𝑆𝑖 ≈ 𝑓𝐼𝐵,𝑆𝑖𝐺𝑒 and 𝑘𝑆𝑖 ≈ 𝑘𝑆𝑖𝐺𝑒 ≈ 5500. It can also be deduced that 𝑓𝐼𝐵,𝑆𝑖1−𝑥𝐺𝑒𝑥 ≈ 1 for 0 ≤𝑥 ≤ 0.2 . Next, we applied 𝑘𝑆𝑖𝐺𝑒 = 5500 and extracted the value of 𝑓𝐼𝑃,𝑆𝑖0.82𝐺𝑒0.18 from Eq. (4.11), fitting it to our experimental data of 𝑅𝐶𝑃,𝑆𝑖0.82𝐺𝑒0.18. The best fit was determined by minimizing the 53 root-mean-square error and it gives 𝑓𝐼𝑃,𝑆𝑖0.82𝐺𝑒0.18 = 0.95 . Referring to Figure 4.9, it is clear that a small vacancy contribution ( 𝑓𝑉𝑃,𝑆𝑖0.82𝐺𝑒0.18 = 0.05 ) is responsible for the U-shape of the 𝑅𝐶𝑃,𝑆𝑖0.82𝐺𝑒0.18. In this work, we aimed to minimize the number of empirical parameters needed to fit the diffusivities in known cases, and, hopefully, bolster our confidence in the results concerning cases yet to be explored experimentally. Figure 4.11 Modeling of carbon effect on P and B diffusivity in Si and strained SiGe. Literature data are from the same references as in Figure 4.8. 4.6 Comparison of carbon impact on different dopants Now that the carbon impact models on defects and dopant diffusivities were established and applied to B and P in Si and SiGe, it is natural to ask whether these models apply to other dopants 00.20.40.60.810 0.05 0.1 0.15 0.2 0.25 0.3 0.35P, 18%Ge (this work)P, Si (Ref.)P, Si (this work)B, 10~20%GeB, Sik = 5500, FI = 1k = 5500, FI = 0.95RCSubstitutional Carbon Concentration (Atom %)𝑹𝑪 =𝟏𝟏 + 𝒌 [𝑪𝒚] 𝑹𝑪 =𝟎. 𝟗𝟓𝟏 + 𝒌 [𝑪𝒚] + 𝟎. 𝟎𝟓 ∙ (𝟏 + 𝒌 [𝑪𝒚]) 54 such as arsenic (As) and antimony (Sb). Although these models are based on data from HBT base-like structures and annealing conditions, they may also work for a wider range of structures and annealing conditions. In this section, we will discuss the applications of the RC model and its implications on As and Sb dopant diffusion. Only one literature work on As and Sb diffusion in Si:C was found, Ref. [48], where the diffusion structures are 150 nm Si/150 nm Si:C/Si. One side of the As and Sb profiles are inside the Si:C, and the diffusion anneals were at 900 ⁰C for 2 to 6 hours. We compared the literature 𝑅𝐶𝐴𝑠,𝑆𝑖 and 𝑅𝐶𝑆𝑏,𝑆𝑖 with the calculations using our RC model. We found that the experimental data agree with our model using 𝑘𝑆𝑖 ≈ 5500, 𝑓𝐼𝐴𝑠,𝑆𝑖 = 0.4 and 𝑓𝐼𝑆𝑏,𝑆𝑖 = 0.02 , the latter two of which are well accepted values derived from numerous literature data [22]. This comparison shows that our models can also apply to other dopants and to non-HBT like structures and annealing conditions. Further experiments can be performed to investigate the dynamics of carbon-defect interactions with varied structure design and time-dependent annealing conditions. Here we converted the carbon molar fraction into relative interstitial concentration ([𝐼]̃/[𝐼]∗), using 𝑘𝑆𝑖 ≈ 5500 and Eq. (4.8), and plotted it against the dopant diffusivity ratio 𝐷/𝐷∗ in a log-log scale so that the trend is more clear (see Figure 4.12). The choice of using [𝐼]̃/[𝐼]∗ instead of using carbon concentration is to include other defect engineering approaches such as thermal nitridation. It shows that this model can effectively catch the carbon impact behavior for all these dopants, i.e. B, P, As, and Sb. The carbon impact on Sb in the original reference [48] was given as 8 ± 2 (not shown in this Figure), which is close to the prediction of our model. Comparing Figure 4.9 and Figure 4.11, one can see that the carbon-introduced diffusivity change is much more sensitive to 𝑓𝐼 when 𝑓𝐼 is close to 1. Considering the equivalent position of interstitial and vacancy in Eq. (2.4), we can describe this feature as “the significant change in diffusivity amplifies the small difference in the non-dominant fractional contribution of diffusion”. Therefore, for vacancy dominant diffusion such as Sb in Si, the suppression of vacancy 55 concentrations can be used to determine 𝑓𝑉𝑆𝑏,𝑆𝑖 in a more accurate way, as 𝑓𝑉𝑆𝑏,𝑆𝑖 is close to 1. One can imagine that if there exists one element that can suppress vacancy concentrations effectively, which is the opposite of carbon, then we would be able to use them as a pair to disturb diffusion and derive the 𝑓𝐼 and 𝑓𝑉 for the diffusing species (similar to the function of oxidation and nitridation which suppress vacancy and interstitial respectively). In this sense, Tin (Sn) is a potential candidate which may serve this purpose. However, this is beyond the scope of our current work and we will not discuss in further detail here. Figure 4.12 Experiment RC data for Sb, As, P and B vs. model predictions. Literature data are from the same references as in Figure 4.8. 0.010.11100.010.11P, 18%GeP, SiB, 11~20%GeB, SiAs, SiSb, SifI = 1.00fI = 0.95fI = 0.40fI = 0.02RC[I] / [I]*[𝐼]̃/[𝐼]∗ 56 4.7 Conclusions In this chapter, the carbon effects on P diffusion in SiGe:C and Si:C under rapid thermal anneal conditions were studied quantitatively. The Ge molar fractions studied are 0 and 18%, and the carbon molar fractions are up to 0.32%. The results showed that the carbon retardation effect is less effective for 𝑆𝑖0.82𝐺𝑒0.18 than Si. In 𝑆𝑖0.82𝐺𝑒0.18: 𝐶, there is an optimum carbon content at around 0.05% to 0.1%, beyond which more carbon incorporation does not add more retardation to P diffusion. These phenomena can be explained by the increased vacancy-mediated diffusion fraction 𝑓𝑉𝑃 in strained SiGe as Ge content increases. Compared to boron, P diffusion in strained Si0.82Ge0.18 has a bigger vacancy-mediated diffusion fraction by about 5%. Empirical models were established to calculate the carbon impact on the time-averaged point defect concentrations and effective diffusivities as a function of carbon content. The strain effect due to carbon incorporation was shown to be negligible. The models were shown to be consistent with experiments on antimony and arsenic diffusion in Si:C and boron diffusion in SiGe:C. These empirical models are easy to calibrate and implement in process simulators. These results also indicate that carbon is very effective in perturbing point defect concentrations, which in turn may provide a useful approach to investigate dopant diffusion and self-diffusion mechanisms. 57 The effect of thermal nitridation on phosphorus diffusion in SiGe and SiGe:C As discussed in 1.1.3, thermal nitridation is known to suppress B and P diffusion in Si, but there are no published records on its industry applications. R&D on SiGe PNP-HBTs provides a good opportunity for the application of this method, as the dopants in PNP-HBTs, B and P, are both mainly interstitial diffusers. It is of technological interest to investigate whether this method is effective. If it is the case, it is a relatively low cost and simple approach to be adopted by industry, as the only change is in the annealing ambient. These questions will be addressed in Chapter 5. 5.1 Literature review 5.1.1 Physical effects of thermal nitridation The term “thermal nitridation” was historically used in several types of silicon wafer processing technologies which involve the chemical element nitrogen [81]–[83]. In this work, we refer to thermal nitridation or nitridation with regard to only “annealing silicon wafers in ammonia (NH3) ambient”. This process involves the reaction of ammonia gas with either 𝑆𝑖𝑂2 (oxynitridation) or the bare silicon surface (direct nitridation). The early studies of thermal nitridation mainly focused on how to produce high-quality thermal silicon nitride films [81], [84], [85]. Hayafuji et al. first investigated the physical phenomena during thermal nitridation of silicon, such as lattice defects [86]. They studied the shrinkage and growth of preexisting oxidation-induced stacking faults (see Figure 5.1) during thermal nitridation of silicon with and without an oxide film. They observed that stacking faults in silicon without an oxide film shrink linearly with nitridation time while stacking faults in silicon with an oxide film grew during nitridation (see Figure 5.2). It was also found that the stacking faults’ shrinkage and growth rates depend on the partial pressure of ammonia. Based on these results, a 58 model was suggested that the shrinkage and growth of stacking faults were caused by the emission and absorption of silicon self-interstitials by the Frank partial dislocations [87]. It also implies that thermal nitridation with silicon cause the change in intrinsic point defect concentrations. This theory was backed by Fahey et al. who studied the impact of thermal nitridation on dopant diffusion [17]. Figure 5.1 An example of the etch figures of the stacking faults taken in a wafer after thermal nitridation at 1100 ℃ for 4 hours [86]. 59 Figure 5.2 The shrinkage of stacking faults during nitridation of silicon without oxide film under an ammonia pressure of 1 kg/cm2 [86]. 5.1.2 Effect of thermal nitridation on dopant diffusion It was found that direct nitridation of the silicon surface retarded boron and phosphorus diffusion while oxynitridation enhanced diffusion of both impurities [17]. The responses of phosphorus, boron, antimony and arsenic diffusion to identical external conditions were compared [88], [18]. The results were interesting in that some impurities exhibit the opposite behavior to other impurities under identical conditions. During direct nitridation where bare silicon wafers were annealed in ammonia, phosphorus and boron diffusion was retarded while antimony diffusion was enhanced. In contrast, for oxynitridation, phosphorus and boron diffusion was enhanced while antimony diffusion was retarded. These results revealed that direct nitridation and oxynitridation cause different changes to the silicon wafers and different impurities have various responses to the same change. Coupled with the theory proposed by Hayafuji, a full picture was drawn that direct nitridation and oxynitridation caused the point defect concentrations to change in opposite ways, 60 and subsequently each impurities’ diffusivity was affected corresponding to their specific diffusion mechanism. Although these arguments could not prove a one-to-one correspondence between macroscopic impurity profiles and microscopic diffusion mechanism, it is still useful for us to understand the various experimental results within a consistent theoretical framework. It also provided convincing clues on the dominant diffusion mechanism of each impurity, which is of great technological importance. A systematic review was given by Fahey et al. with a focus on the relation between intrinsic point defects and impurity diffusion [35]. The experimental results were explained by assuming Sb diffusion to be dominated by the vacancy mechanism, P and B diffusion to be dominated by the interstitial mechanism, and As to have comparable components of both types of mechanisms. The early studies in the 1980s were conducted mainly using spreading resistance analysis to obtain impurity profiles and using analytic arguments to derive results. Ural et al. offered improvements by implementing the SIMS technique and performing numerical simulations [36]. The results of impurity diffusion mechanism were given with six sets of numerical simulations with various levels of approximations and assumptions on point defect dynamics. These results generally agreed with previous studies and confirmed the hypotheses on defect injection conditions and impurity diffusion mechanisms. Another issue of thermal nitridation experiments is the kinetics of injected point defects. Mogi et al. reported a time-dependent study using doped superlattices and demonstrated the transport capacity of injected point defects [19]. They used antimony diffusion profiles as the indicator to reveal the vacancy supersaturation during thermal nitridation as a function of time and temperature. They showed that vacancies diffused through the entire epilayer (650nm) within 60 minutes at 860℃ or 15 minutes at 910℃. 61 5.1.3 Thermal nitridation on Si and SiGe Several studies have been reported on thermal nitridation effects on both Si and SiGe. Karunaratne et al. studied the effect of thermal nitridation on boron diffusion in Si, strained Si0.89Ge0.11 and Si0.89Ge0.11:C and the diffusion suppression factors in the different materials were compared [16] [89]. The results showed that the extent of B diffusion suppression in Si and Si0.89Ge0.11 is comparable between 940 ℃ and 1050 ℃. The nitridation-induced suppression factor in the carbon-incorporated sample was weaker than in a carbon-free sample. This indicated that the effects of thermal nitridation and carbon incorporation can both suppress boron diffusion and their impacts can be added to some extent, although the total effect is not simply the sum of the two effects. Bonar et al. studied the effect of thermal nitridation on Sb diffusion in Si and SiGe [90]. They showed that the Sb diffusion enhancement factor was smaller in Si0.9Ge0.1 than in Si. Given the present understanding of point defect injection conditions, it is interesting to explore its technological potential to control dopant diffusion. R&D on SiGe PNP-HBTs provides a good opportunity for the application of this method, as the dopants in PNP-HBTs, B and P, are both interstitial diffusers, and may be retarded at the same time with the vacancy injection method. However, there is no experimental study available on whether this method works for P in strained-SiGe systems and to what extent. It is also not clear whether this effect will work with other approaches (e.g. carbon incorporation) and to what extent. If this method is effective, it is a relatively low cost and simple approach that can be adopted by industry, as the only change is in the annealing ambient. These questions are addressed in this work. 5.2 Diffusion experiments We designed P in SiGe and SiGe:C structures and annealed them in three different defect injection conditions, including the thermal nitridation (vacancy injection) condition. The effectiveness of the thermal nitridation method was studied based on these experiments. Apart from 62 the industry applications, thermal nitridation also provides an important approach to reveal P diffusion mechanisms in SiGe and SiGe:C. 5.2.1 Sample structures and annealing conditions design To understand the thermal nitridation effect on P diffusion behavior and mechanisms in 𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦, a wafer matrix with various Ge and C contents was designed and grown epitaxially. Figure 5.3 shows the epitaxial structure design, where 𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦 are sandwiched between two Si layers. To be relevant to PNP SiGe HBT applications, the x values were chosen to be 0.09 and 0.18 and the y values ranged from 0 to 0.09%. The naming convention of the structures is to use the Ge and C percentages to name a 𝑆𝑖/𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦/𝑆𝑖 structure. Table 5.1 Epitaxial wafer matrix used in this work. Structure Wafer name Ge mole fraction Carbon mole fraction Si0.91Ge0.09 Ge09C0 9% 0 Si0.91Ge0.09:C Ge09C009 9% 0.09% Si0.82Ge0.18 Ge18C0 18% 0 Si0.82Ge0.18:C Ge18C006 18% 0.06% Si0.82Ge0.18:C Ge18C009 18% 0.09% As discussed in 3.2.2, all the wafers were grown on 8 inch (100) Czochralski p-type Si wafers in an ASM Epsilon 2000 RPCVD reactor. Ammonia anneals were performed using an AccuThermo AW610 atmospheric rapid thermal processing (RTP) system. The target annealing 63 temperatures was 980 to 1000⁰C, and the annealing time ranged from 35 seconds to 6 minutes. The thermal annealing conditions were designed such that P diffusion happened inside the fully strained-SiGe or fully strained-SiGe:C layers. Figure 5.3 Schematic diagram of the structure design. 5.2.2 Defect injection using masking layers To extract the thermal nitridation effect, control samples were needed to show diffusion under inert conditions (no defect injection). To avoid the run-to-run temperature variation, we needed to anneal the control samples together with those under thermal nitridation conditions in one RTA anneal. Therefore, for the control samples, we used SiO2/SiNx as masking layers, which has been proven to protect the underlying structure and create equivalently inert diffusion conditions, even when ammonia is present [67]. An oxide layer was also used without nitride to create interstitial injection conditions due to the oxynitridation reaction in ammonia. The oxide and nitride masking layers were deposited by plasma-enhanced chemical vapour deposition (PECVD). Figure 5.4 (a) and (b) show the cross sections of bare samples and masked control samples. By annealing the 64 three types of samples in the same RTA run, we avoided run-to-run variation, and made the samples directly comparable to show the effect of vacancy injection and interstitial injection. The annealing conditions were designed carefully to avoid the strain relaxation in SiGe and SiGe:C layers. After annealing, the masking layers were etched away by 1:10 hydrogen fluoride (HF) containing buffered oxide etch (BOE) solution, and SIMS was used to obtain P, Ge and C concentration profiles. Figure 5.4 Structures used for different defect injection conditions during thermal nitridation experiments. 5.3 Experimental results 5.3.1 Strain and composition analysis (113) reciprocal space mapping was performed on as-grown and annealed samples at room temperature. Figure 5.5 shows the RSMs of some annealed samples, which were annealed with the largest thermal budget for the corresponding structure. These results clearly show that the in-plane lattice constants of the 𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦 layer are the same as that of the Si substrate after annealing. This means that these two layers are fully strained to the substrate during diffusion, and no strain relaxation was observed for the structures after diffusion. These results confirm that 65 P diffusion happened in fully strained SiGe or fully strained SiGe:C. In the following discussion, SiGe stands for fully strained SiGe unless otherwise noted. 𝜔 − 2𝜃 rocking scans of (004) Bragg reflection were performed to confirm the substitutional carbon fraction of each structure (Figure 5.6). The measurement procedures are the same as described in Ch.4.3.2. Figure 5.5 (113) reciprocal space mapping of annealed samples with the largest thermal budget for each type of structures: a) Ge09C0; b) Ge09C009; c) Ge18C0; d) Ge18C009. 66 Figure 5.6 XRD (004) scan of SiGe and SiGe:C samples used for thermal nitridation experiments. The shifts of (004) peaks show Ge and carbon molar fraction differences. 5.3.2 Concentration profiling and P diffusivity ratio extraction Figure 5.7 shows the SIMS profiles of as-grown and annealed samples for each structure. Both the as-grown and annealed P peaks are within the constant Ge and C region. The shape of the diffused P peak is close to a Gaussian distribution, indicating that the concentration dependence of P diffusivity in SiGe is not strong. C diffusion is seen to be much faster than P and Ge diffusion. Figure 5.7(a) shows that the P diffusion under the nitridation condition was retarded compared to the inert condition, and it was enhanced under the oxynitridation condition as expected, which showed that the defect engineering methods by capping layers were valid. As Si-Ge interdiffusion is much slower than P diffusion, annealed Ge profiles are very close to the as-grown Ge profiles; 67 therefore, only Ge as-grown profiles are shown for clarity. In sample Ge09C009, shown in Figure 5.7(b), there is no clear difference between P diffusion under the nitridation condition and that under the inert condition. However, carbon diffusion is retarded under the nitridation condition, as seen in the carbon rising edge at the depth below 40 nm and in the falling edge in the 140 to 180 nm depth range. As carbon is an interstitial diffuser in Si, this is a clear indication of the interstitial undersaturation caused by vacancy injection. The carbon profiles in all other diffusion structures are similar, and are not shown for clarity. It is also shown that this interstitial undersaturation penetrates to a depth of at least 200 nm. Figure 5.7(c) to (e) shows the results from samples with 18% Ge but different C contents. Interestingly, in Ge18C006 and Ge18C009, the P profile under the nitridation condition is not showing any retardation, but rather a slight enhancement in P diffusion is seen. 68 69 Figure 5.7 SIMS profiles of as-grown and annealed samples and TSUPREM-4 fittings. The solid lines on top of the P profiles are the TSUPREM-4 fittings. The nitridation impact factor on P diffusivity, 𝑅𝑁𝑇𝐷 =𝐷𝑃𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝐷𝑃𝑖𝑛𝑒𝑟𝑡 , was extracted by fitting SIMS profiles and listed in Table 5.2. As the temperature accuracy of the RTP system is estimated to be a few tens of degrees, some temperature adjustments were needed to fit the SIMS profiles. These adjustments do not affect the nitridation impact factor, as each set of samples for comparison 70 were annealed in the same run. The error bar of 𝑅𝑁𝑇𝐷 =𝐷𝑃𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝐷𝑃𝑖𝑛𝑒𝑟𝑡 is estimated to be ±20% due to the SIMS error and TSUPREM-4 fitting errors. Figure 5.7 includes some examples of the best fitting curves and the corresponding SIMS profiles for 𝑅𝑁𝑇𝐷 extractions. The nitridation impact factor, 𝑅𝑁𝑇𝐷 =𝐷𝑃𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝐷𝑃𝑖𝑛𝑒𝑟𝑡 , are shown in Figure 5.8 as a function of carbon molar fraction. A value greater than one indicates a diffusivity enhancement, whereas smaller than one means retardation. In reality, during annealing, the point defect concentrations, and therefore the diffusivity ratios, are functions of time. Therefore, the diffusivity values extracted in this work are time-averaged values. From Figure 5.8, we can see that thermal nitridation is still effective in retarding P diffusion in SiGe for up to 18% Ge, but when 0.06% or 0.09% carbon is present, this retardation is no longer effective. As a comparison, relevant literature data are also plotted in this Figure (Ref.[16] and [36]). Similar to the P case, 𝑅𝑁𝑇𝐷 for B (𝑅𝑁𝑇𝐷,𝐵) increases with the C concentration for Si and Si0.89Ge0.11. However, 𝑅𝑁𝑇𝐷,𝐵 is less than 1 with C, which means that even with C, the nitridation still retards B diffusion on top of the carbon retardation. While for Si0.82Ge0.18:C samples, 𝑅𝑁𝑇𝐷,𝑃 can be larger than 1, showing that the thermal nitriation of Si increases P diffusion. 71 Table 5.2 Nitridation impact factor for phosphorus RNTD,P under non-equilibrium point defect conditions caused by nitridation. The annealing temperatures quoted here are the temperature extracted using TSUPREM-4 fitting to account for uncertainty. Wafer name Annealing condition* Nitridation impact factor 𝑅𝑁𝑇𝐷,𝑃 Ge09C0 973 ⁰C 1 min 0.45 ± 0.09 Ge09C009 920 ⁰C 6 min 1.00 ± 0.20 Ge18C0 982 ⁰C 35 s 0.79 ± 0.16 Ge18C006 982 ⁰C 35 s 1.41 ± 0.28 Ge18C009 952 ⁰C 60 s 1.19 ± 0.24 0.20.40.60.811.21.41.61.80 0.02 0.04 0.06 0.08 0.1RNTDCarbon Concentration (Atom %)P in 18%GeP in 9%GeB in 11%Ge (Ref.)B in Si (Ref.)P in Si (Ref.) Figure 5.8 Nitridation impact factor RNTD as a function of carbon % at various Ge %. Error bars of the data without carbon are not shown for clarity. The literature data quoted in this figure are obtained from Ref.[36] for P in Si and from Ref.[16] for B in Si and SiGe. 72 5.4 Discussions 5.4.1 Major observations Two major correlations can be observed from the experimental results shown in Figure 5.8. First, for a certain Ge molar fraction, 𝑅𝑁𝑇𝐷 goes up with the increasing C content, i.e., the effectiveness of P diffusivity retardation under the vacancy injection condition decreases with the increasing carbon content. For the carbon molar fractions we used, 0.06% and 0.09%, no further retardation effect is seen from the thermal nitridation. Second, for a certain C molar fraction, 𝑅𝑁𝑇𝐷 goes up with the increasing Ge fraction. These effects can be explained using the point defect assisted diffusion theory. To analyze the thermal nitridation effect on P diffusion we can apply the point-defect assisted diffusion theory, as described in Ch.2. In this chapter, we define the impact of thermal nitridation as 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒 ≡𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒:𝐶 = 𝑓𝐼[𝐼] [𝐼]∗+ 𝑓𝑉[𝑉] [𝑉]∗ . (5.1) In SiGe:C samples, 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒:𝐶 shows only the thermal nitridation effect, not the overall carbon plus thermal nitridation effect. 5.4.2 Ge fraction dependence of thermal nitridation impact factor On the Ge fraction dependence, the observation in Figure 5.8 is consistent with the results shown in Ch.4. It was shown that when the matrix changes from Si to Si0.82Ge0.18, the vacancy mediated diffusion fraction 𝑓𝑉 increases from 0 to 0.05. Therefore, in Si, 𝑓𝑉𝑆𝑖 = 0, we have 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖 =[𝐼] [𝐼]∗. (5.2) In Si0.82Ge0.18, 𝑓𝑉𝑆𝑖0.82𝐺𝑒0.18 = 0.05, and 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖0.82𝐺𝑒0.18 = 0.95[𝐼] [𝐼]∗+ 0.05[𝑉] [𝑉]∗. (5.3) The experiment results in Figure 5.8 shows that 73 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖0.82𝐺𝑒0.18 > 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖091𝐺𝑒0.09 > 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖 . (5.4) Although we do not have direct evidence to show that [𝐼] [𝐼]∗ values in Si0.82Ge0.18 and Si0.91Ge0.09 are close to that in Si under the same thermal nitridation condition, 𝑓𝑉𝑆𝑖0.82𝐺𝑒0.18 > 𝑓𝑉𝑆𝑖0.91𝐺𝑒0.09 > 𝑓𝑉𝑆𝑖, and [𝐼] [𝐼]∗Si0.82Ge0.18≈[𝐼] [𝐼]∗Si0.91Ge0.09≈[𝐼] [𝐼]∗Si is one possible explanation of the Ge dependence observed. When carbon is present, the carbon concentration is low and introduces interstitial undersaturation. We assume that carbon concentration studied does not change the P diffusion mechanisms, i.e., 𝑓𝑉𝑆𝑖0.82𝐺𝑒0.18 = 𝑓𝑉𝑆𝑖0.82𝐺𝑒0.18:𝐶, and 𝑓𝑉𝑆𝑖0.91𝐺𝑒0.09 = 𝑓𝑉𝑆𝑖0.91𝐺𝑒0.09:𝐶. Therefore, it is expected that the Ge dependence doesn’t change with the carbon addition, as seen in Figure 5.8: 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖0.82𝐺𝑒0.18:C > 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖091𝐺𝑒0.09:C (5.5) 5.4.3 Carbon fraction dependence of thermal nitridation impact factor On the carbon fraction dependence, the results indicate that when both are present, the nitridation retardation effect is weaker, i.e, 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒:𝐶 > 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒 . For P in Si0.82Ge0.18:C, enhancement in P diffusion is seen instead of retardation, which may be a result of high vacancy supersaturation acting on the non-negligible 𝑓𝑉 diffusion term. Therefore, the method of using thermal nitridation to control B and P diffusion only works for Si and very low C and Ge content SiGe:C, which is the area below 𝑅𝑁𝑇𝐷 = 1 line in Figure 5.8. The trend of the increased 𝑅𝑁𝑇𝐷,𝑃 with the carbon presence agrees with the 𝑅𝑁𝑇𝐷,𝐵 results from Ref.[16]. However, 𝑅𝑁𝑇𝐷,𝐵 < 1 still holds, where 𝑅𝑁𝑇𝐷,𝑃 can be larger than 1. This difference may be a result of 𝑓𝑉,𝑃𝑆𝑖𝐺𝑒 > 𝑓𝑉,𝐵𝑆𝑖𝐺𝑒 , as suggested from Ref.[91] and the above Ge dependence discussion. As generally accepted, both the carbon incorporation and the thermal nitridation of Si act on dopant diffusion by reducing the concentrations of interstitials [11], [12], [69], [86]. Historically, 74 the influence of thermal nitridation on point defects was observed indirectly by the growth or shrinkage of stacking faults and dopant diffusion [86], as the point defect concentrations are too low to be measured directly. When both methods are applied simultaneously, the nitridation impact is less, i.e., 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒:𝐶 > 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒 . As P is not a pure interstitial diffuser in SiGe for the Ge fractions studied, while B is, using B diffusion as an indirect measurement of the point defect concentration change when both carbon method and nitridation method are combined, although outside the scope of this work, will be valuable to give a clearer defect concentration picture. 5.4.4 Thermal nitridation and equivalent [𝑰] [𝑰]∗ in Si0.82Ge0.18 In Ch.4, we have systematically studied the carbon impact on defect engineering in Si and SiGe, and expressed the carbon impact factor 𝑅𝐶,𝑃𝑆𝑖𝐺𝑒:𝐶 as a function of the carbon content and the interstitial undersaturation ratio [𝐼] [𝐼]∗, which is less than 1. A natural thought is to compare 𝑅𝐶,𝑃 and 𝑅𝑁𝑇𝐷,𝑃, and relate the latter with an equivalent carbon content and [𝐼] [𝐼]∗. The relation between the interstitial undersaturation ratio [𝐼]∗ [𝐼] in Si:C and Si0.82Ge0.18:C and the carbon molar fraction [C] was established by (equivalent as Eq.(4.8)) [𝐼]∗ [𝐼]= 1 + 5500 [𝐶] (5.6) By definition: 𝑅𝐶,𝑃𝑆𝑖𝐺𝑒:𝐶 ≡𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 . (5.7) Also, the steady-state condition holds, 𝑅𝐶,𝑃𝑆𝑖0.82𝐺𝑒0.18:𝐶 = 0.95[𝐼] [𝐼]∗+ 0.05[𝑉] [𝑉]∗= 0.95[𝐼] [𝐼]∗+ 0.05 [𝐼]∗ [𝐼] (5.8) By definition, 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒:𝐶 ≡𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒:𝐶 =𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 ∗𝑅𝐶,𝑃 (5.9) 75 Therefore, 𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 = 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖𝐺𝑒:𝐶 ∗ 𝑅𝐶,𝑃𝑆𝑖𝐺𝑒:𝐶 (5.10) 𝑅𝐶,𝑃𝑆𝑖𝐺𝑒:𝐶 has been measured in Si0.82Ge0.18:C in Ref. [91], and the annealing in that work was performed in nitrogen ambient at 1000 ⁰C for 15 to 120 seconds, close to the conditions in this work. Therefore, it is reasonable to assume 𝑅𝐶,𝑃𝑆𝑖𝐺𝑒:𝐶 obtained from Ref. [91] can be applied to this work. Combining (11) and (13), we can calculate [𝐼] [𝐼]∗ for SiGe:C under thermal nitridation. P diffusivity under a defect engineering condition divided by the equilibrium P diffusivity, i.e., 𝐷𝑃,𝑑𝑒𝑓𝑒𝑐𝑡−𝑒𝑛𝑔𝑖𝑛𝑒𝑒𝑟𝑒𝑑𝐷𝑃∗ is plotted as a function of the interstitial undersaturation ratio [𝐼] [𝐼]∗ in Si0.82Ge0.18 in Figure 5.9. Depending on which defect engineering method(s) is (are) used, 𝐷𝑃,𝑑𝑒𝑓𝑒𝑐𝑡−𝑒𝑛𝑔𝑖𝑛𝑒𝑒𝑟𝑒𝑑𝐷𝑃∗ can be 𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 or 𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 or 𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 . Table 5.3 Summary of the impact factor of carbon, thermal nitridation, and both methods on P diffusion in Si0.82Ge0.18 and the corresponding interstitial undersaturation ratio [I]/[I]*. The impact of carbon is adopted from Ref. [91]. Carbon Concentration 𝑅𝐶,𝑃𝑆𝑖0.82𝐺𝑒0.18:𝐶 [𝐼] [𝐼]∗ (carbon effect only) 𝑅𝑁𝑇𝐷,𝑃𝑆𝑖0.82𝐺𝑒0.18:𝐶 𝐷𝑛𝑖𝑡𝑟𝑖𝑑𝑎𝑡𝑖𝑜𝑛𝑆𝑖𝐺𝑒:𝐶𝐷𝑖𝑛𝑒𝑟𝑡𝑆𝑖𝐺𝑒 for 𝑆𝑖0.82𝐺𝑒0.18 [𝐼] [𝐼]∗ (carbon + nitridation effect) 0 1 1 0.79 ± 0.16 0.79 ± 0.16 0.77 0.06% 0.44 0.200 1.41 ± 0.28 0.62± 0.12 0.094 0.09% 0.46 0.165 1.19 ± 0.24 0.55± 0.11 0.114 From Table 5.3 we can see that thermal nitridation further injected vacancies and suppressed self-interstitial concentrations by 31% to 53%. These results are close to the observation of vacancy 76 injection impact on boron diffusion in strained SiGe:C with 11% Ge from Ref.[16], which provided extra 8% to 42% interstitial suppression on top of the carbon impact for the temperature range of 940 to 1000 ⁰C. 00.511.500.20.40.60.8100.511.5DP,defect-engineered/D* PInterstitial Undersaturation [I]/[I]*fI = 0.95, fV = 0.05 Figure 5.9 P diffusivity under a defect engineering condition divided by the P diffusivity ratio as a function of the interstitial undersaturation ratio [I]/[I]* in Si0.82Ge0.18. It is based on local equilibrium assumption of [I][V] = [I]*[V]*.This figure was adapted from Ref. [91]. To further investigate the effects of thermal nitridation on point defects and dopant diffusion in strained SiGe, although beyond the scope of this work, one may consider to use the comprehensive procedure previously used for the study of dopant diffusion in Si [36]. This procedure combines the nitridation and oxidation processes, and uses identical conditions for multiple dopants, such as boron and antimony, where a parallel analysis of different dopants can be performed in a more rigorously way to reveal key parameters such as [𝐼] [𝐼]∗ and [𝑉] [𝑉]∗. 77 5.5 Conclusions To conclude, we reported an experimental study of the thermal nitridation effects on phosphorus (P) diffusion in strained 𝑆𝑖1−𝑥𝐺𝑒𝑥 and strained 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝐶𝑦 with up to 18 at.% Ge and 0.09 at.% carbon. P diffusivities under thermal nitridation (vacancy injection) and effective inert conditions were compared. The results show that thermal nitridation can retard P diffusion in SiGe, but the effectiveness of this retardation method decreases with increasing Ge and C content. When 0.06% or 0.09% carbon is present in Si0.82Ge0.18, thermal nitridation slightly increases P diffusivity compared to the inert conditions. The Ge dependence can be explained by the increasing contribution from the vacancy-assisted mechanism for P diffusion in strained SiGe with the increasing Ge content. In terms of interstitial undersaturation ratio [𝐼] [𝐼]∗ in Si0.82Ge0.18, thermal nitridation can further decrease [𝐼] [𝐼]∗ by 31 to 53% on top of the carbon effect. 78 Segregation of phosphorus in SiGe After the discussion of P diffusion in uniform Si and SiGe with a focus on diffusivities and diffusion mechanisms, let us turn to another issue: the dopant distribution across graded SiGe layers. Although SiGe layers with a constant Ge molar fraction are essential for diffusion studies, graded SiGe layers are more commonly used in SiGe HBTs. This involves another universal phenomenon regarding the redistribution of a chemical element in an inhomogeneous solid material system besides diffusion, which is called segregation. In strained SiGe systems, phosphorus tends to segregate from higher Ge regions into lower Ge regions, while boron has exactly the opposite tendency. This is especially a problem for thin SiGe layers in PNP SiGe HBTs because segregation causes broadened P profiles, which contradicts the goal to achieve high and narrow P profiles for high-speed PNP HBTs. Segregation phenomenon cannot be described by a simple modification in the diffusivity using Fick’s law alone. As the segregation happens simultaneously with diffusion, we have to account for the segregation-related transport together with the diffusion transport. Without knowledge and models of segregation, one can not predict dopant profiles accurately. 6.1 Introduction Let us clarify some terms before further discussions. In a broader sense, diffusion refers to the net mass transport in a material system. Many phenomena can happen with diffusion such as grain growth, phase change, segregation, etc. To align our studies with the conventions of literature on diffusion in semiconductors, we restrict the term “diffusion” to the net mass transport of an impurity in homogeneous solids where the material composition can be considered to be uniform. In contrast, the term “segregation” is used to describe the impurity distribution in inhomogeneous solids where the material composition is not constant. The 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝑃 material systems provide good examples of the “coupled diffusion and segregation” phenomenon. By using a colon in 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝑃, we 79 indicate that only the concentrations of Si and Ge affect the material property of SiGe, as a binary alloy, and impurity P has negligible impact except for providing charge carriers as a dopant. That being said, only when Ge molar fraction “x” changes, do we refer to this material as inhomogeneous. In many practical structures, segregation and diffusion are decoupled in space, where segregation happens and is modeled only at interfaces. For Si and SiGe material systems, most dopant segregation studies have been limited to sharp interfaces between two different materials, such as Si/SiO2 and Si/Si3N4. Segregation is typically modeled using a flux balance equation with a segregation coefficient in many commercial simulation tools such as TSUPREM-4TM. It should be noted that when there is not a sharp interface, but a gradual change in the solid-state solvent, such as a gradual change in Ge concentration in a SiGe alloy (see Figure 6.1), segregation also happens, and it is coupled with dopant diffusion. In the Sentaurus ProcessTM [78], the segregation model for a graded solid solvent is available, but the calibration was based only on two sets of experiments using sharp SiGe interfaces under equilibrium conditions, which are not relevant to typical industry structures or annealing conditions. Figure 6.1 Illustration of sharp and graded interfaces. SiGe HBTs are good examples to investigate coupled diffusion and segregation. The epitaxial strained heterostructures provide ideal testbeds. A typical SiGe HBT base layer has an epitaxial SiGe layer with a graded Ge slope on one or both sides, and the layer is doped with either boron (B) for a NPN type HBT or with phosphorus (P) for a PNP type HBT. During fabrication steps with 80 high temperatures, dopant diffusion and segregation happen simultaneously in the same space. Segregation is beneficial for NPN HBTs as B prefers to stay in SiGe base layers with higher Ge content, while in PNP HBTs, P segregates out of SiGe base layers to Si layers and makes the vertical doping profile hard to control. In this case, coupled dopant diffusion and segregation models are required to describe the dopant profile evolution across the Ge concentration slopes. There are only a handful of experimental studies on dopant segregation in SiGe, which were typically performed on sharp Si/SiGe interfaces using long annealing times [15], [74], [92]. Due to the low spatial resolution of the concentration profiles, these data can hardly be used for model calibration purposes. Therefore, more accurate and industry relevant experimental data for segregation model calibration are in great need. In this study, we re-derived a diffusion-segregation model clearly showing the diffusion and segregation fluxes, and performed P segregation experiments at SiGe slopes to verify the model. Compared to segregation experiments using sharp SiGe/Si interfaces, where one segregation coefficient can be measured for one SiGe composition, segregation experiments using graded SiGe provide a new approach to calibrate and/or extract the segregation coefficients for a range of SiGe compositions, and the experimental structures are directly relevant to the SiGe HBT industry. Our coupled diffusion-segregation model was then calibrated using these, which enables more accurate prediction of dopant profile evolution during diffusion. The theoretical derivations of coupled diffusion and segregation model of this work is generic to all inhomogenous solid solutions, and the experiments and model calibration are specifically relevant to PNP SiGe HBTs. 6.2 Literature review Early studies on dopant segregation in SiGe material were mainly focused on how to introduce certain amounts of doping into elemental semiconductors for the development of fabrication technology [93]–[95]. This research focused mainly on the segregation coefficients (called “distribution coefficients” at that time) for various situations together with chemical reactions 81 during crystal growth, often from a thermodynamic equilibrium perspective [93], [96] [97]. Segregation, as a part of mass transport phenomena during high temperature annealing after the structure fabrication, has been investigated along with the development of SiGe epitaxy techniques. Hu et al. performed segregation experiments where dopants diffused from a polycrystalline silicon layer on top of the Si/SiGe/Si epitaxial layers [98]. It was found that boron tends to segregate into the SiGe layer while phosphorus and arsenic tend to segregate away from it (see figure 6.2 as an example). More experimental studies were reported with drive-in diffusion ([92], [99]) and in-situ doping ([15], [100]) conditions. Figure 6.2 Comparison of phosphorus concentration profiles in Ge-doped (solid curve) and control (dotted curve) silicon. Figure reprinted from Ref. [98] with permission. From a theoretical perspective, Hu laid the groundwork by analyzing the chemical potential of a dopant in inhomogeneous media [101][102]. In an inhomogeneous material, the driving force for the diffusion of a dopant is no longer simply given by its concentration gradient, but by the gradient of its chemical potential. The chemical potential of a dopant, considering the atomic system and electronic subsystem, is given by [102]: 82 𝜇1 = 𝑘𝑇𝑙𝑛 (𝑁12𝑁𝐿(𝑁𝐶𝑁𝑉)1/2) − 𝑍1Δ𝑊𝑖(𝑁2) − Δ𝐸𝑏1(𝑁2) +12𝐸𝑔(𝑁2) + 𝜃𝛽1𝛽2𝑁2 + 𝜃𝛽12𝑁1 (6.1) where 𝑁1, 𝑁2 and 𝑁𝐿 are the concentrations of the dopant, Ge and lattice sites, respectively. The bonding energy Δ𝐸𝑏1 , the intrinsic work function 𝑊𝑖 , the bandgap 𝐸𝑔 , the lattice contraction/expansion coefficient 𝛽2, the conduction- and valence-band density of states 𝑁𝐶 and 𝑁𝑉 are all functions of the germanium concentration 𝑁2. Theoretical calculations suggested that B segregation is dominated by stress effects, while P and As is dominated by changes in electronic band structure [103]. With approximation of low doping concentration for 𝑁1 < 1 × 1020𝑐𝑚−3, the dopant flux 𝐽1 was given by the gradient of chemical potential as: 𝐽1 = −𝐷1𝑁1𝑘𝑇𝜕𝜇1𝜕𝑥 = −𝐷1𝜕𝑁1𝜕𝑥+𝐷1𝑁12[1𝑁𝐶𝜕𝑁𝐶𝜕𝑁2+1𝑁𝑉𝜕𝑁𝑉𝜕𝑁2+2𝑘𝑇(𝑍1𝜕𝑊𝑖𝜕𝑁2+𝜕𝐸𝑏1𝜕𝑁2−12𝜕𝐸𝑔𝜕𝑁2− 𝜃𝛽1𝛽2)]𝜕𝑁2𝜕𝑥 . (6.2) With the gradients 𝜕𝑁1𝜕𝑥 and 𝜕𝑁2𝜕𝑥, the diffusion and segregation flux can be calculated simultaneously. However, this model requires many parameters, such as the effective density of states (𝑁𝐶 and 𝑁𝑉), lattice contraction coefficients (𝛽1 and 𝛽2) and several local energy terms, which are difficult to obtain from theory or experiments. On the other hand, You et al. derived a simpler equation based on thermodynamic principles (parameter labels are modified here for consistency) [104] : 𝐽1 = − 𝐷1 (𝜕𝑁1𝜕𝑥−𝑁1𝑁1𝑒𝑞𝜕𝑁1𝑒𝑞𝜕𝑥) , (6.3) where 𝑁1𝑒𝑞 was defined as “normalized thermal equilibrium concentration”. It was given by [105]: 𝑁1𝑒𝑞𝑁𝐿= exp (−𝑔𝑓𝑘𝐵𝑇) , (6.4) 83 where 𝑁𝐿 was lattice site concentration and 𝑔𝑓 was the formation Gibbs free energy of the diffusing species. However, there is no clear explanation on the precise meaning of 𝑁1𝑒𝑞. For example, 𝑁1𝑒𝑞 was explained as solubility and was used to define segregation coefficient 𝑘𝑠𝑒𝑔 as [106]: 𝑘𝑠𝑒𝑔 ≡𝐶𝛽A1𝐶𝛽A2 =𝐶𝛽A1𝑠𝑜𝑙𝑢𝑏𝑖𝑙𝑖𝑡𝑦𝐶𝛽A2𝑠𝑜𝑙𝑢𝑏𝑖𝑙𝑖𝑡𝑦 . (6.5) By using solid solubility ratios, this definition implies that the segregation coefficient is a constant, independent of dopant concentration, which is theoretically inappropriate [102]. On the other hand, the definition of “solubility” is also ambiguous since it has at least two different meanings: 1) as described by Henry’s law: the solubility of a gas in a liquid is directly proportional to the partial pressure of the gas above the liquid; 2) the maximum possible concentration of an element which can be dissolved into a solid state material without creating second phase. Therefore, the parameters’ definitions were confusing and it was not clear on how to use Eq. (6.3) for modeling of coupled diffusion-segregation phenomenon. From these references, one can see that a gap between experimental results and theoretical modeling still exists. In the following sections, we will address this issue through detailed derivation and definitions and specially designed experimental verifications. 6.3 Derivation of coupled of diffusion-segregation equation In some classic diffusion textbooks, such as the one by Shewmon [107], we can find some discussions on ternary alloys, where the mass flux of one element is not only driven by its own concentration gradient 𝜕𝑁1𝜕𝑥, the normal diffusion term, but is also driven by the gradient of the second element 𝜕𝑁2𝜕𝑥, which is the segregation term. SiGe:P can be considered as a ternary alloy. Let’s at first review the derivations of the mass fluxes, which are based on classic thermodynamics 84 and previous work [23]. We use a more convenient definition of 𝑘𝑠𝑒𝑔, and separate the diffusion flux and segregation flux in the data fitting and simulations in Section 6.4. 6.3.1 The driving force for diffusion in solid state materials In solid-state materials, diffusion and segregation are statistical results of random jumps of atoms. This phenomenon can be described thermodynamically as if the atoms are driven by the chemical potential gradient of the diffusing species. The mass flux can be expressed as 𝐽 = 𝑁 ∙ ?⃗? = 𝑁 ∙ (𝑀 ∙ 𝐹 ) (6.6) where N, V⃗ , 𝑀 and F⃗ denote concentration, velocity, mobility and the driving force, respectively [108]. For a one dimension case along the z axis, the driving force is opposite to the chemical potential gradient, 𝐹 ∙ 𝑒𝑧 = −𝜕𝜇𝜕𝑧 , (6.7) where 𝜇 is the chemical potential and 𝑒𝑧 is the unit vector of z axis. As pointed out by Atkins [109], such a force should be considered as a “thermodynamic force”. It is not necessarily a real force pushing the particles down the slope of the chemical potential. Considering the fact that the solid state diffusion is microscopically the random jump of atoms, this force may represent the spontaneous tendency of the atoms to disperse as a consequence of the Second Law of thermodynamics and the hunt for maximum entropy. Usually, this force is so called the driving force for diffusion. 6.3.2 Derivation of coupled diffusion-segregation equation Now let us consider a three-element system, e.g. 𝑆𝑖1−𝑥𝐺𝑒𝑥: 𝑃. Let 𝑁1 be the concentration of P, the solute, and 𝑁2 be the concentration of Ge, which represent the material composition. Here, 𝑁1 is small enough so that its impact on material property is negligible and 𝑁2 is large enough to cause material property change. To differentiate their different roles, we put a colon “:” between 85 𝑆𝑖1−𝑥𝐺𝑒𝑥 and 𝑃. For generality purpose, we will treat P as a solute 𝛽 and 𝑆𝑖1−𝑥𝐺𝑒𝑥 as the solution A, respectively. Then the chemical potential gradient of solute 𝛽 in solution A can be expressed as 𝜕𝜇𝛽A𝜕𝑧=𝜕𝜇𝛽A𝜕𝑁1𝜕𝑁1𝜕𝑧+𝜕𝜇𝛽A𝜕𝑁2𝜕𝑁2𝜕𝑧 , (6.8) where the first term on the right-hand side is the traditional diffusion driving force term as in the Fick’s law, and the second term is the driving force of segregation. As discussed in Ch.6.1, in this thesis, “diffusion” is used in its narrower sense, which means the net transport of substance due to its chemical potential gradient in a given homogeneous diffusion media. Segregation is to describe the net transport due to inhomogeneous diffusion media. By separating the segregation term, users of a simulation tool can evaluate the magnitude of the segregation term to see if it is significant enough to be included for the price of longer computation time. In the case of P in biaxially-fully-strained SiGe, we treat 𝑁2 (𝑁2 = 𝑥𝐺𝑒 ∙ 𝑁𝑡𝑜𝑡𝑎𝑙) as the sole representative of the material property. It can be approximated that 𝜕𝑁𝑡𝑜𝑡𝑎𝑙𝜕𝑁2≈ 0 when Ge concentration is low enough so that its impact on lattice concentration is negligible. This approximation generally holds in this study. Combining Eq.(6.6)~(6.8), we obtain 𝐽 = −𝑁1 ∙ 𝑀 ∙ (𝜕𝜇𝛽A𝜕𝑁1𝜕𝑁1𝜕𝑧+𝜕𝜇𝛽A𝜕𝑁2𝜕𝑁2𝜕𝑧) . (6.9) In a homogeneous structure, the material is uniform, i.e. 𝜕𝑁2𝜕𝑧= 0, and there is no segregation. When 𝜕𝑁2𝜕𝑧≠ 0 , we need to find out the expression of 𝜕𝜇𝛽A𝜕𝑁2 , which represents the non-uniform nature of inhomogeneous structure. Based on classic thermodynamics, at a fixed pressure and temperature, the chemical potential of solute 𝛽 in solution A, such as P in P doped SiGe, is 𝜇𝛽A = 𝜇𝛽A,⦵ + 𝑅𝑇𝑙𝑛(𝑎𝛽A), and 𝑎𝛽A = 𝑥𝛽A𝛾𝛽A , 𝑥𝛽A =𝑁1𝑁𝑡𝑜𝑡𝑎𝑙, (6.10) 86 where 𝜇𝛽A,⦵, R, T , 𝑎𝛽A , 𝑥𝛽A, 𝛾𝛽A, 𝑁1 and 𝑁𝑡𝑜𝑡𝑎𝑙 denote the standard chemical potential of 𝛽 in solid state material A, ideal gas constant, the absolute temperature, the activity, the molar fraction, the chemical activity coefficient of 𝛽 in A, the number of 𝛽 atoms, and the total number of atoms in solution A. The chemical potential of 𝛽 with respect to the Ge concentration can now be expressed as 𝜕𝜇𝛽A𝜕𝑁2= 𝑅𝑇𝜕𝑙𝑛(𝑥𝛽A∙𝛾𝛽A(𝑁2))𝜕𝑁2 . (6.11) To find out the expression of 𝜕𝜇𝛽A𝜕𝑁2, let us consider the situation when solute 𝛽 is distributed within two adjacent and different solutions (or phases) A1 and A2. The concept of thermal equilibrium in a system requires that the chemical potential of each species involved is constant in the system, i.e. equal in adjacent phases. This condition can be expressed as 𝜇𝛽A1(𝑥𝛽A1) = 𝜇𝛽A2(𝑥𝛽A2) . (6.12) Using the relation of Eq.(6.10), it can be further expressed as 𝜇𝛽A1,⦵ + 𝑅𝑇𝑙𝑛(𝛾𝛽A1) + 𝑅𝑇𝑙𝑛(𝑥𝛽A1) = 𝜇𝛽A2,⦵ + 𝑅𝑇𝑙𝑛(𝛾𝛽A2) + 𝑅𝑇𝑙𝑛(𝑥𝛽A2) . (6.13) Then the segregation coefficient 𝑘𝑠𝑒𝑔 is defined as 𝑘𝑠𝑒𝑔 (𝑁2𝑁𝑡𝑜𝑡𝑎𝑙) ≡𝑥𝛽A2𝑥𝛽A1=𝛾𝛽A1𝛾𝛽A2∙ exp (𝜇𝛽A1,⦵ −𝜇𝛽A2,⦵ 𝑅𝑇) . (6.14) Since the Si, SiGe and Ge discussed here are all with the same crystal structure, the standard chemical potential of 𝛽 in these materials are identical [23] (𝜇𝛽Si,⦵ = 𝜇𝛽SiGe,⦵ = 𝜇𝛽Ge,⦵ ) . Therefore the segregation coefficient corresponds to the inverse ratio of activity coefficient 𝑘𝑠𝑒𝑔 ≡𝑁1A2𝑁1A1 =𝛾𝛽A1𝛾𝛽A2 . (6.15) Also, the concentration of 𝛽 can be considered as independent of 𝑁2 (i.e. 𝜕𝑙𝑛(𝑥𝛽A)𝜕𝑁2= 0). Then Eq. (6.11) can be expressed as 𝜕𝜇𝛽A𝜕𝑁2= 𝑅𝑇𝜕𝑙𝑛(𝛾𝛽A(𝑁2))𝜕𝑁2= −𝑅𝑇𝜕𝑙𝑛(𝑘𝑠𝑒𝑔(𝑁2))𝜕𝑁2 . (6.16) 87 Combining Eq. (6.8), (6.10) and (6.16), the chemical potential gradient of 𝛽 can be expressed as 𝜕𝜇𝛽A𝜕𝑧=𝑅𝑇𝑁1∙𝜕𝑁1𝜕𝑧− 𝑅𝑇𝜕 ln(𝑘𝑠𝑒𝑔)𝜕𝑁2𝜕𝑁2𝜕𝑧 . (6.17) Applying this equation into Eq.(6.9), we get the mass flux expression as: 𝐽 = −𝑁1 ∙ 𝑀 [𝑅𝑇𝑁1∙𝜕𝑁1𝜕𝑧− 𝑅𝑇𝜕 ln(𝑘𝑠𝑒𝑔)𝜕𝑁2𝜕𝑁2𝜕𝑧] = −𝐷1𝜕𝑁1𝜕𝑧+ 𝐷1𝑁1𝜕 ln(𝑘𝑠𝑒𝑔)𝜕𝑁2𝜕𝑁2𝜕𝑧 . (6.18) On the right-hand side, the first term is the traditional diffusion flux term as in the Fick’s law, and the second term is the segregation flux. With a chosen reference material, if 𝑘𝑠𝑒𝑔 changes with N2 (such as Ge fraction in SiGe alloys), then segregation happens in an inhomogeneous media with N2 being non-uniform. Mathematically, Eq. (6.18) can be written as 𝐽 = −𝐷1 ∙1𝑘𝑠𝑒𝑔∙𝜕(𝑁1𝑘𝑠𝑒𝑔)𝜕𝑧 , (6.19) where both diffusion and segregation terms are integrated into one spatial derivative. This is the format of segregation model suggested in Sentaurus Process Advanced Calibration Kit, where the segregation effect was included in the “stress model”, as stress also causes segregation. It is, however, useful for a simulation tool to have these two terms separated, as in Eq. (6.18). In many process simulation tools, diffusion flux is well modelled, where segregation flux is not commonly included except at material interfaces. By separating these two terms, one can do a rough estimation on the relative magnitude of both terms and decide whether the segregation flux is significant enough to be included, at the cost of longer simulation time. During the above derivations, P in SiGe is used as an example, and Eq. (6.18) is generic to all coupled diffusion and segregation phenomena in alloys due to alloy composition change only. Stress can also change the chemical potential of a species and thus drive segregation, which should add another term in Eq. (6.18). As many semiconductor alloy systems are pseudomorphic, stress is a function of the alloy composition. Therefore, the stress effect can be included in kseg. In other words, if kseg in Eq. (6.18) is the segregation coefficient for fully strained material, i.e. 𝑘𝑠𝑒𝑔(𝑁2) =𝑁1A2,fully strained𝑁1A1 , then Eq. (6.18) can be used for pseudomorphic alloy systems. In SiGe HBTs, SiGe 88 layers are designed to be fully strained, where Eq. (6.18) applies. In this work, the segregation coefficients studied are for Si/fully-strained-SiGe. 6.4 Experimental verification 6.4.1 Sample structure design and experiments Now that the model was derived, experiments were designed to validate our model. To observe the coupled diffusion segregation phenomena, we designed structures to let dopants diffuse and segregate across the graded SiGe regions. The peak concentrations of Ge used were 14% and 18% in this work (see Figure 6.3). The P peak was designed to be inside the triangle SiGe region so that we could observe P segregation out from the Ge triangle to Si [98]. Figure 6.3 Schematic diagram of the structure design. The Ge peak concentrations used were 14% and 18%. As discussed in Ch.3.2.2, all the wafers were grown on 8 inch (100) Czochralski p-type Si wafers in an ASM Epsilon 2000 RPCVD reactor. Inert anneals at 900 ⁰C were performed in nitrogen ambient using an enclosed Linkam TS1200 high temperature heating stage. The thermal budgets of the annealing were selected to avoid strain relaxation. At 900 ⁰C, the ramp-up rate was 89 about 60 ⁰C/minute and ramp-down rate was about 150 ⁰C/minute. The temperature ramp rates were fast enough such that we can neglect the thermal budget in the ramp-up and ramp-down stages. The temperature uncertainty of the heating stage is less than ± 10 ⁰C, which is about ± 50% in terms of P diffusivity at the annealing temperature. This uncertainty was taken into account in our simulations. Depth profiles of P and Ge were measured by SIMS for as-grown and annealed samples. Ion yield and sputter rate were calibrated to accommodate the matrix composition change. 6.4.2 Parameters for modeling As discussed in Chapter 3, our first choice was to use TSUPREM-4 for modeling and simulations. As it does not allow the incorporation of the 2nd gradient term in its flux calculations, we needed to write our own numerical calculation code. The computer coding for the new algorithm was done in Matlab. To implement the model according to Eq. (6.18), a finite difference approach was used to numerically calculate P diffusion and segregation fluxes simultaneously and the P profile evolution with time. For this purpose, we need to determine two parameters, the P diffusivity 𝐷1 and the segregation coefficient as a function of Ge concentration, 𝑘𝑠𝑒𝑔(𝑁2). Since these two parameters are both temperature dependent, we chose a fixed temperature of 900 ⁰C for our annealing experiments as a starting point. The P diffusivity 𝐷1 is based on the calibrated values for Si including the Fermi effect [22]. The Ge fraction (𝑥𝐺𝑒) dependence of 𝐷1 was also taken into account according to the data from Ref.4. For P diffusivity in SiGe up to 18% Ge, P diffusivity increases monotonically with the Ge fraction, and the difference between P in Si and P in Si0.82Ge0.18 is about 40%, which is not big. Therefore, we approximated the P diffusivity in SiGe with a simple linear function of Ge molar fraction. A small factor 𝐷𝑠𝑐𝑎𝑙𝑒 is applied to accommodate the experimental temperature uncertainty. Then the P diffusivity 𝐷1 is expressed as: 𝐷1 = 𝐷𝑠𝑐𝑎𝑙𝑒 × (1 + 𝑥𝐺𝑒) × [𝐷0 + 𝐷−(𝑛𝑛𝑖) + 𝐷=(𝑛𝑛𝑖)2] , (6.20) 90 where 𝑥𝐺𝑒 = 𝑁2/𝑁𝑡𝑜𝑡𝑎𝑙 . This expression includes the concentration-dependent diffusivity, as described in Chapter 2. The parameters such as 𝐷0 , 𝐷−and 𝐷= are a set of experimentally-calibrated valued adopted by TSUPREM-4 (see Table 6.1). 𝑛𝑖 is the intrinsic carrier concentration in Si and can be expressed by 𝑛𝑖 = 3.87 × 1016 × exp (−0.605𝑒𝑉𝑘𝐵𝑇) cm-3. The value of 𝑛 is calculated according to local charge neutrality condition. Table 6.1 List of major parameters used in TSUPREM-4 simulations. The diffusivities are in μm2/min. Concentrations are in cm-3. Parameter Prefactor Activation energy (eV) 𝑫𝟎 2.31 × 1010 3.66 𝑫− 2.66 × 1010 4.0 𝑫= 2.652 × 1011 4.37 𝒏𝒊 3.87 × 1016 0.605 To determine the values of 𝑘𝑠𝑒𝑔 as a function of 𝑥𝐺𝑒, we used the available literature data as a starting point. Figure 6.4 is a summary of experimental data where Ge concentration is expressed in at. %. It shows that 𝑘𝑠𝑒𝑔(𝑥𝐺𝑒) is not sensitive to temperature within the range of 800⁰C to 950⁰C at the low Ge end. This may be due to two possible reasons: 1) the temperature dependence is insignificant compared to the experimental accuracy; and 2) segregation coefficients were extracted as a ratio of P concentration on SiGe side over that on Si side from “on-going” diffusion-segregation profiles in these data. From the time evolution of P profiles (as shown later in this article), we can see that this ratio is not a constant. This data should be considered showing a “trend of the segregation phenomenon” rather than the determination of the value of 𝑘𝑠𝑒𝑔 . Meanwhile, the 91 experiment condition of some experiments may not be consistent, e.g., the thickness of SiGe layers may be too thick to be dislocation-free [92]. In principle, the P segregation coefficient can be expressed by an Arrhenius relation where the activation energy is a function of Ge content [101], [110]. 𝑘𝑠𝑒𝑔 = exp (−𝑥𝐺𝑒∙𝐸𝑠𝑒𝑔𝑘𝑇) (6.21) Here, 𝐸𝑠𝑒𝑔 the activation energy, which indicates how strong the Ge affects P segregation. In the following simulations, the values of 𝐸𝑠𝑒𝑔 were determined by fitting to the experimental data. 00.20.40.60.810 0.2 0.4 0.6 0.8 1950oC (Hu)900oC (Christensen)900oC (Rucker)800oC (Kobayashi)ksegGe Molar Fraction Figure 6.4 Experimental data of P segregation coefficient in SiGe. The data points in diamonds are from Ref. [92], which were obtained from SiGe films thicker than 0.2 micron where no strain analysis was available. The remaining data are extracted from Ref. [15], [74], [98]. 6.4.3 SIMS data and data fitting The SIMS profiles of four different structures and/or annealing conditions and the corresponding fitting results are plotted in Figure 6.3 (a) ~ (d). As-grown SIMS profiles of P and 92 Ge were taken as pre-annealed profiles and post-annealing profiles were numerically calculated. As mentioned, any annealing temperature uncertainty will influence the diffusivity, which was treated with a fitting parameter 𝐷𝑠𝑐𝑎𝑙𝑒. The fitting parameters, 𝐷𝑠𝑐𝑎𝑙𝑒 and 𝐸𝑠𝑒𝑔, were determined when the best match to the post-annealing SIMS profiles was achieved. For each diffused profile, 𝐷𝑠𝑐𝑎𝑙𝑒 determines the amount of diffusion, seen as the broadening of the annealed profiles, while 𝐸𝑠𝑒𝑔 determines the shape of the diffused peak, i.e., the dip inside the SiGe and the peak in the surrounding Si layers. These two effects are independent, so the best fitting 𝐷𝑠𝑐𝑎𝑙𝑒 and 𝐸𝑠𝑒𝑔 can be obtained for each diffused profile. Figure 6.3 (a) to (d) compare the as-grown and annealed Ge and P SIMS with the best-fitting calculation for the four structure and condition combinations. It can be seen that our model can accurately describe the coupled diffusion-segregation behavior of P across the gradient 𝑆𝑖1−𝑥𝐺𝑒𝑥 region. In comparison, the diffusion-only profiles cannot catch the segregation phenomena. The best-fitting 𝐸𝑠𝑒𝑔 and the corresponding segregation coefficient (𝑘𝑠𝑒𝑔) for each set of data and segregation coefficients from literature data are summarized in Figure 6.4. From this comparison, within the temperature range of 800 ⁰C to 950 ⁰C and experimental errors, we can see that 𝐸𝑠𝑒𝑔 =0.5 𝑒𝑉 is a good compromise for all the available studies. 93 10171018101910200510150 50 100 150900 oC 20min14% Ge peakSIMS_As-grownSIMS_AnnealedDiffusion-segregationDiffusion onlyGermaniumP Concentration (cm-3)Ge Molar Fraction (%)Depth (nm)(a)Eseg = 0.55 eV10171018101910200510150 50 100 150900 oC 10min14% Ge peakSIMS_As-grownSIMS_AnnealedDiffusion-segregationDiffusion onlyGermaniumP Concentration (cm-3)Ge Molar Fraction (%)Depth (nm)(b)Eseg = 0.29 eV 94 Figure 6.5 Ge and P SIMS profiles and best-fitting curves using the coupled diffusion-segregation model in Eq. (6.18). Simulations with diffusion-only model were plotted for comparison. 1017101810191020051015200 20 40 60 80 100 120 140 160900 oC 20min18% Ge peakSIMS_As-grownSIMS_AnnealedDiffusion-segregationDiffusion onlyGermaniumP Concentration (cm-3)Ge Molar Fraction (%)Depth (nm)(c)Eseg = 0.56 eV1017101810191020051015200 40 80 120 160900 oC 10min18% Ge peakSIMS_As-grownSIMS_AnnealedDiffusion-segregationDiffusion onlyGermaniumP Concentration (cm-3)Ge Molar Fraction (%)Depth (nm)(d)Eseg = 0.49 eV 95 00.20.40.60.810 0.05 0.1 0.15 0.2 0.25 0.3950oC (Hu)900oC (Christensen)900oC (Rucker)800oC (Kobayashi)Eseg = 0.29 eVEseg = 0.49 eVEseg = 0.55 eVEseg = 0.56 eVksegGe Molar Fraction Figure 6.6 The best-fitting Eseg and the corresponding segregation coefficients extracted from coupled diffusion-segregation experiments. The dashed parts are extrapolation using the Eseg values. Literature data are from the same reference as in Figure 6.4. 6.4.4 Model prediction and P profile time evolution With the model and the calibrated parameter Eseg, we can simulate the time evolution of the P profile during a segregation-diffusion process, as seen in Figure 6.7. The pre-annealed P and Ge profiles are the as-grown Ge and P profiles in the 14% Ge structure. We can see that the segregation effect is more and more obvious as time increases, which is a dynamic process. It should be noted that in previous segregation coefficient studies [74], [92], the ratio of P concentration in SiGe over that in Si was measured and quoted as the segregation coefficient. As we can see from Figure 6.5 (b) and (c), this ratio also depends on time. Under the assumption that SiGe layer stays fully strained 96 during the annealing, it takes time for P segregation to reach a stable value, and only after that can one use the concentration ratio as the segregation value. The difficulty in doing that for Si/SiGe system is that when the time is too long, SiGe layers may start to relax or interdiffuse with Si, and the segregation is then influenced by the stress relaxation and interdiffusion. Therefore, it makes more sense to use the coupled segregation and diffusion model to extract the segregation coefficients from experiments with annealing time short enough to avoid stress relaxation. The above experiments and analysis demonstrated a new approach in measuring segregation coefficients. As the coupled diffusion segregation model calculates dopant fluxes, it is capable to calculate dopant profile evolution during short anneals before the concentration ratio reaches a stable value. The segregation coefficient extraction is done by fitting the entire dopant profiles instead of by calculating the concentration ratios at a sharp interface. The use of this model eliminates the need of long time annealing and thick layers to extract segregation coefficients, which are not practical for many heterostructures with lattice mismatch strains. Future studies are necessary to fine-tune our model for P diffusion in SiGe and to apply our method to different material systems. 97 101710181019102010210510150 50 100 150SIMS As-grownSIMS Annealed4 min8 min12 min16 min20 minGeP Concentration (cm-3) Ge Molar Fraction (%)Depth (nm)(a)10161017101810191020051015200 200 400 600P initialP 10 hrsGe initial Ge 10 hrsP Concentration (cm-3) Ge Molar Fraction (%)Depth (nm)900 oCGe(b) 98 Figure 6.7 (a) Simulated time evolution of P profiles during annealing in comparison with SIMS data. The peak Ge fraction is 14%, (b) simulated time evolution of P profiles across sharp Si/SiGe interfaces during annealing and (c) extracted P concentration ratio at either side of the left Si/SiGe interface as a function of time (from 1 to 10 hours). The initial P and Ge have box-like profiles. The kseg for Si/Si0.8Ge0.2 (fully strained) is shown for comparison. 6.5 Impact of P segregation on HBT device performance At this point, one may wonder if the segregation effect is significant enough to impact the device performance of a PNP HBT. To address this question, we performed one-dimensional device simulations with a mainstream device simulation tool MediciTM [111] using their default models for PNP HBTs. The initial P, C and Ge profile design was modified from the NPN HBT design in Ref.[5]. We chose a 3 nm thick box-like P as-grown profile and calculated the annealed profiles using Eq. (6.18). Two annealing conditions (5 s and 10 s at 1000 ⁰C) were simulated. Two types of P base doping profiles were loaded in Medici for comparison: one was simulated with the 00.20.40.60.810 2 4 6 8 10P Concentration RatioTime (hours)ksegCP,Si / CP,Si0.8Ge0.2900 oC(c) 99 segregation-diffusion model and the other was simulated with the diffusion model only, as shown in Figure 6.8. The electrical characteristics were simulated using industry-proven drift-diffusion model. For HBT device performance, depending on the applications, two figures of merits are important. For HBTs in low-speed applications such as the low-speed input stages to amplifier cells, the common-emitter current gain (𝛽) is very important; while the common-emitter short-circuit cut-off frequency 𝑓𝑇 is less important. In radio-frequency (RF) applications, 𝑓𝑇 is a very important figure of merit. It is desired to achieve the fastest transit times possible in HBTs, therefore they are operated at large current densities to minimize impacts of junction capacitances. The impact of the segregation effect on 𝛽 of this PNP HBT is shown in Figure 6.9, where the collector voltage was set to 2 V. It shows that the segregation effect reduces the current gain over a wide range of collector current. The impact of segregation effect on 𝑓𝑇 is shown in Figure 6.10. It showed that in the high collector current IC window, where the RF applications usually operate, the segregation effect reduces the cut-off frequency 𝑓𝑇 greatly. These results show that the segregation effect has a strong impact on PNP HBT performance. It is therefore important to use the coupled diffusion-segregation model in the base doping and Ge profile design. 100 10161017101810191020051015202520 30 40 50 60 70 801000 oC5s w/o segregation5s w/ segregation10s w/o segregation10s w/ segregationGeNet Doping Concentration (cm-3)Ge Molar Fraction (%)Depth (nm)Emitter(P)Base(N)Collector(P) Figure 6.8 Doping and Ge profiles used for one-dimensional HBT device simulation. Figure 6.9 One-dimensional device simulation results showing the impact of the P segregation effect on PNP HBT current gain (𝜷). The P base profiles are simulated using the diffusion-segregation model and the diffusion only model for comparison. The annealing conditions used are 5 s and 10 s at 1000 ⁰C. 110100-12.0 -10.0 -8.0 -6.0 -4.0Current Gain (β)Log (IC) (A/um)5s w/o segregation5s w/ Segregation10s w/o segregation10s w/ Segregation 101 Figure 6.10 One-dimensional device simulation results showing the impact of the P segregation effect on PNP HBT cut-off frequency (fT). The P base profiles are simulated using the diffusion-segregation model and the diffusion only model for comparison. The annealing conditions used are 5 s and 10 s at 1000 ⁰C. 6.6 Conclusions In this chapter, a coupled diffusion and segregation model was re-derived, where the contributions from diffusion and segregation to dopant flux are explicitly shown. The model is generic to coupled diffusion and segregation in inhomogeneous alloys, and provides a new approach in segregation coefficient extraction. This model is especially helpful for heterostructures with lattice mismatch strains. Experiments of coupled P diffusion and segregation were performed with graded SiGe layers for Ge molar fractions up to 0.18, which are relevant to PNP SiGe HBTs. The model was shown to catch both diffusion and segregation behavior well. The diffusion-segregation model for P in SiGe alloys was calibrated and 𝐸𝑠𝑒𝑔 = 0.5 𝑒𝑉 is suggested for the temperature range from 800 to 950 ⁰C. 891011-6.0 -5.0 -4.0 -3.0Log (fT ) (Hz)Log (IC) (A/um)5s w/o segregation5s w segregation10s w/o Segregation10s w/ Segregation 102 The above experiments and analysis demonstrated a new approach in measuring segregation coefficients, which eliminates the needs of long time annealing and thick layers, which are not practical for many heterostructures with lattice mismatch strains. 103 Contributions and future work 7.1 Contributions of this work To our best knowledge, this work is the first systematic and quantitative experimental investigation of phosphorus diffusion and segregation phenomena in strained SiGe and SiGe:C systems. The findings are relevant to device applications with SiGe, SiGe:C and Si:C materials such as SiGe HBTs. 1) The experimental data represents the first systematic investigation of the dependence of P diffusion in SiGe:C due to C incorporation. The effective range and limitations of using carbon to control P diffusion in SiGe were shown, which was not studied previously in literature. The germanium dependence of the carbon impact was also studied. This finding has important technological impact since it implies that there is a design trade-off between Ge concentration and P profile sharpness. These results are of great interest for the development of PNP SiGe HBTs with steep P profiles. 2) Scientifically, we demonstrated an approach to study the diffusion mechanisms 𝑓𝐼 and 𝑓𝑉 with carbon induced defect engineering, which has not used before. We developed a point defect-based semi-empirical model to explain the P diffusion behavior at various carbon and Ge concentrations. It was found that the there is an increasing contribution of vacancy-assisted mechanism toward P diffusion with increasing Ge content. The interstitial diffusion fraction 𝑓𝐼 of P was shown to be 95% in Si0.82Ge0.18. These models are much easier to calibrate and implement than the dopant-defect reaction and transport equations. 3) To our best knowledge, this work is also the first experimental study investigating thermal nitridation effect on P diffusion in SiGe and SiGe:C. The result shows that thermal nitridation can retard P diffusion in SiGe with up to 18% Ge content, but the effectiveness of this retardation decreases with increasing Ge and C content. These phenomena can be explained by the increasing 104 contribution from the vacancy-assisted mechanism for P diffusion in strained SiGe with increasing Ge content. 4) We developed a coupled diffusion-segregation model to simulate the segregation phenomenon in inhomogeneous SiGe media. The definition of segregation coefficient was clarified. By experimental verification, we demonstrated that this model is capable of simulating the P segregation phenomenon across gradient SiGe regions. This model is generic to coupled diffusion and segregation in inhomogeneous alloys, and provides a new approach in segregation coefficient extraction, which is especially helpful for heterostructures with lattice mismatch strains. Since the segregation phenomenon causes undesired P profiles for PNP HBT device, this model is of practical importance for the device design trade-off consideration for optimized performance. 7.2 Suggestions for future work Rather than completing, this work opens even more interesting research topics which may be explored in the future. 1) The understanding of the native point defects in silicon has been pursued for about eight decades, while the effort for this mission is still ongoing [12], [84]. For SiGe materials, the challenge is even greater with an extra dimension of complexity. For example, the identity of self-interstitial (could be Si or Ge) becomes a question and their roles in dopant diffusion is still not clear. Strained SiGe is inherently not a standalone bulk material, thus imposing greater difficulty in experimentation. Further studies may be performed, either first-principle calculations or new experimental methodology, to understand atomic diffusion mechanisms in details. 2) This work studied the segregation across SiGe layers without carbon, and the annealing was performed with a heating stage. The natural extension of this work is to study the segregation across SiGe:C layers with an industry standard RTA tool. Due to the constraint of our heating stage, we were not able to reach 1000 to 1100 ⁰C range, which is typically used in the industry. 105 3) The segregation model is generic to any inhomogeneous material system with continuous matrix change. It can be further applied to other dopants and host materials, like in III-V semiconductors. It will be worthwhile to study other material systems and compare the results to SiGe systems. 106 Bibliography [1] D. J. Paul, “Si/SiGe heterostructures: from material and physics to devices and circuits,” Semiconductor Science and Technology, vol. 19, no. 10, pp. R75–R108, Oct. 2004. [2] P. Ashburn, SiGe Heterojunction Bipolar Transistors. John Wiley & Sons, 2003. [3] Y. Shi and G. Niu, “Vertical profile design and transit time analysis of nano-scale SiGe HBTs for Terahertz fT,” in Bipolar/BiCMOS Circuits and Technology, 2004. Proceedings of the 2004 Meeting, 2004, pp. 213 – 216. [4] Y. Shi and G. Niu, “2-D analysis of device parasitics for 800/1000 GHz fT/fmax SiGe HBT,” in Bipolar/BiCMOS Circuits and Technology Meeting, 2005. Proceedings of the, 2005, pp. 252 – 255. [5] D. Bolze, B. Heinemann, J. Gelpey, S. McCoy, and W. Lerch, “Millisecond annealing of high-performance SiGe HBTs,” in 17th International Conference on Advanced Thermal Processing of Semiconductors, 2009. RTP ’09, 2009, pp. 1 –11. [6] J. D. Cressler, “Issues and Opportunities for Complementary SiGe HBT Technology,” ECS Trans., vol. 3, no. 7, pp. 893–911, 2006. [7] D. M. Monticelli, “The future of complementary bipolar,” in Bipolar/BiCMOS Circuits and Technology, 2004. Proceedings of the 2004 Meeting, 2004, pp. 21 – 25. [8] D. Knoll, B. Heinemann, K. E. Ehwald, A. Fox, H. Rucker, R. Barth, D. Bolze, T. Grabolla, U. Haak, J. Drews, B. Kuck, S. Marschmeyer, H. H. Richter, M. Chaimanee, O. Fursenko, P. Schley, B. Tillack, K. Kopke, Y. Yamamoto, H. E. Wulf, and D. Wolansky, “A Low-Cost, High-Performance, High-Voltage Complementary BiCMOS Process,” in Electron Devices Meeting, 2006. IEDM ’06. International, 2006, pp. 1 –4. [9] B. Tillack, B. Heinemann, D. Knoll, H. Rücker, and Y. Yamamoto, “Base doping and dopant profile control of SiGe npn and pnp HBTs,” Appl. Surf. Sci., vol. 254, no. 19, pp. 6013–6016, Jul. 2008. [10] “The International Technology Roadmap for Semiconductors.” 2013. [11] R. Scholz, U. Gösele, J.-Y. Huh, and T. Y. Tan, “Carbon-induced undersaturation of silicon self-interstitials,” Applied Physics Letters, vol. 72, no. 2, pp. 200–202, Jan. 1998. [12] H. Rücker, B. Heinemann, W. Röpke, R. Kurps, D. Krüger, G. Lippert, and H. J. Osten, “Suppressed diffusion of boron and carbon in carbon-rich silicon,” Appl. Phys. Lett., vol. 73, no. 12, pp. 1682–1684, Sep. 1998. [13] H. Rucker, B. Heinemann, D. Bolze, D. Knoll, D. Kruger, R. Kurps, H. J. Osten, P. Schley, B. Tillack, and P. Zaumseil, “Dopant diffusion in C-doped Si and SiGe: physical model and experimental verification,” in IEEE Int. Electron Devices Meeting Technical Digest pp 345, 1999, pp. 345 –348. [14] H. Rücker and B. Heinemann, “Tailoring dopant diffusion for advanced SiGe:C heterojunction bipolar transistors,” Solid-State Electron., vol. 44, no. 5, pp. 783–789, May 2000. [15] H. Rucker, B. Heinemann, R. Kurps, and Y. Yamamoto, “Dopant Diffusion in SiGeC Alloys,” ECS Trans., vol. 3, no. 7, pp. 1069–1075, Oct. 2006. [16] M. S. A. Karunaratne, A. F. W. Willoughby, J. M. Bonar, J. Zhang, and P. Ashburn, “Effect of point defect injection on diffusion of boron in silicon and silicon–germanium in the presence of carbon,” J. Appl. Phys., vol. 97, no. 11, pp. 113531–113531–7, Jun. 2005. [17] P. Fahey, R. W. Dutton, and M. Moslehi, “Effect of thermal nitridation processes on boron and phosphorus diffusion in 〈100〉 silicon,” Applied Physics Letters, vol. 43, no. 7, pp. 683–685, Oct. 1983. 107 [18] P. Fahey, G. Barbuscia, M. Moslehi, and R. W. Dutton, “Kinetics of thermal nitridation processes in the study of dopant diffusion mechanisms in silicon,” Applied Physics Letters, vol. 46, no. 8, pp. 784–786, Apr. 1985. [19] T. K. Mogi, M. O. Thompson, H.-J. Gossmann, J. M. Poate, and H. S. Luftman, “Thermal nitridation enhanced diffusion of Sb and Si(100) doping superlattices,” Applied Physics Letters, vol. 69, no. 9, pp. 1273–1275, Aug. 1996. [20] A. Ural, “Atomic-scale mechanisms of dopant and self-diffusion in silicon,” Ph.D., Stanford University, United States -- California, 2001. [21] N. E. B. Cowern, S. Simdyankin, C. Ahn, N. S. Bennett, J. P. Goss, J.-M. Hartmann, A. Pakfar, S. Hamm, J. Valentin, E. Napolitani, D. De Salvador, E. Bruno, and S. Mirabella, “Extended Point Defects in Crystalline Materials: Ge and Si,” Phys. Rev. Lett., vol. 110, no. 15, p. 155501, Apr. 2013. [22] J. D. Plummer, M. D. Deal, and P. B. Griffin, Silicon VLSI Technology: Fundamentals, Practice and Modeling. Prentice Hall, 2000. [23] P. Pichler, Intrinsic Point Defects, Impurities, and Their Diffusion in Silicon. Springer, 2004. [24] A. Einstein, B. Podolsky, and N. Rosen, “Can Quantum-Mechanical Description of Physical Reality Be Considered Complete?,” Phys. Rev., vol. 47, no. 10, pp. 777–780, May 1935. [25] B. Pullman, The Atom in the History of Human Thought. Oxford University Press, 2001. [26] H. Mehrer, Diffusion in Solids : Fundamentals, Methods, Materials, Diffusion-Controlled Processes. Berlin/Heidelberg, DEU: Springer, 2007. [27] A. Fick, “V. On liquid diffusion,” Philosophical Magazine Series 4, vol. 10, no. 63, pp. 30–39, Jul. 1855. [28] V. V. Levitin, Atom Vibrations in Solids: Amplitudes and Frequencies. Cambridge Scientific Publishers Limited, 2004. [29] H. B. Huntington and F. Seitz, “Mechanism for Self-Diffusion in Metallic Copper,” Phys. Rev., vol. 61, no. 5–6, pp. 315–325, Mar. 1942. [30] H. B. Huntington and F. Seitz, “Energy for Diffusion by Direct Interchange,” Phys. Rev., vol. 76, no. 11, pp. 1728–1728, Dec. 1949. [31] K. C. Pandey, “Diffusion without Vacancies or Interstitials: A New Concerted Exchange Mechanism,” Phys. Rev. Lett., vol. 57, no. 18, pp. 2287–2290, Nov. 1986. [32] H. Bracht, E. E. Haller, and R. Clark-Phelps, “Silicon Self-Diffusion in Isotope Heterostructures,” Phys. Rev. Lett., vol. 81, no. 2, pp. 393–396, Jul. 1998. [33] H. Bracht, H. H. Silvestri, I. D. Sharp, and E. E. Haller, “Self- and foreign-atom diffusion in semiconductor isotope heterostructures. II. Experimental results for silicon,” Phys. Rev. B, vol. 75, no. 3, p. 035211, Jan. 2007. [34] A. Chroneos and H. Bracht, “Diffusion of n-type dopants in germanium,” Applied Physics Reviews, vol. 1, no. 1, p. 011301, Mar. 2014. [35] P. M. Fahey, P. B. Griffin, and J. D. Plummer, “Point defects and dopant diffusion in silicon,” Rev. Mod. Phys., vol. 61, no. 2, pp. 289–384, Apr. 1989. [36] A. Ural, P. B. Griffin, and J. D. Plummer, “Fractional contributions of microscopic diffusion mechanisms for common dopants and self-diffusion in silicon,” J. Appl. Phys., vol. 85, no. 9, pp. 6440–6446, May 1999. [37] A. Ural, P. B. Griffin, and J. D. Plummer, “Self-Diffusion in Silicon: Similarity between the Properties of Native Point Defects,” Phys. Rev. Lett., vol. 83, no. 17, pp. 3454–3457, Oct. 1999. [38] H. Bracht and E. E. Haller, “Comment on ‘Self-Diffusion in Silicon: Similarity between the Properties of Native Point Defects,’” Phys. Rev. Lett., vol. 85, no. 22, pp. 4835–4835, Nov. 2000. [39] A. Ural, P. B. Griffin, and J. D. Plummer, “Ural, Griffin, and Plummer Reply:,” Phys. Rev. Lett., vol. 85, no. 22, pp. 4836–4836, Nov. 2000. [40] TSUPREM-4TM, version 2007.12, (Synopsys, Inc., Mountain View, CA, 2007), . 108 [41] B. Tillack and Y. Yamamoto, “Atomic Layer Doping for Future Si Based Devices,” ECS Transactions, vol. 22, no. 1, pp. 121–131, 2009. [42] H. Bracht, J. F. Pedersen, N. Zangenberg, A. N. Larsen, E. E. Haller, G. Lulli, and M. Posselt, “Radiation Enhanced Silicon Self-Diffusion and the Silicon Vacancy at High Temperatures,” Phys. Rev. Lett., vol. 91, no. 24, p. 245502, Dec. 2003. [43] H. . Osten, D. Knoll, and H. Rücker, “Dopant diffusion control by adding carbon into Si and SiGe: principles and device application,” Materials Science and Engineering: B, vol. 87, no. 3, pp. 262–270, Dec. 2001. [44] A. Fischer, H.-J. Osten, and H. Richter, “An equilibrium model for buried SiGe strained layers,” Solid-State Electronics, vol. 44, no. 5, pp. 869–873, May 2000. [45] J. W. Matthews and A. E. Blakeslee, “Defects in epitaxial multilayers: I. Misfit dislocations,” Journal of Crystal Growth, vol. 27, pp. 118–125, Dec. 1974. [46] A. Zhylik, A. Benediktovich, A. Ulyanenkov, H. Guerault, M. Myronov, A. Dobbie, D. R. Leadley, and T. Ulyanenkova, “High-resolution x-ray diffraction investigation of relaxation and dislocations in SiGe layers grown on (001), (011), and (111) Si substrates,” Journal of Applied Physics, vol. 109, no. 12, p. 123714, Jun. 2011. [47] H. Nitta, J. Tanabe, M. Sakuraba, and J. Murota, “Carbon effect on strain compensation in Si1−x−yGexCy films epitaxially grown on Si(100),” Thin Solid Films, vol. 508, no. 1–2, pp. 140–142, Jun. 2006. [48] H. Rücker, B. Heinemann, and R. Kurps, “Nonequilibrium point defects and dopant diffusion in carbon-rich silicon,” Phys. Rev. B, vol. 64, no. 7, p. 073202, Jul. 2001. [49] J. P. Liu and H. J. Osten, “Substitutional carbon incorporation during Si1−x-yGexCy growth on Si(100) by molecular-beam epitaxy: Dependence on germanium and carbon,” Applied Physics Letters, vol. 76, no. 24, pp. 3546–3548, Jun. 2000. [50] H. J. Osten, M. Kim, K. Pressel, and P. Zaumseil, “Substitutional versus interstitial carbon incorporation during pseudomorphic growth of Si1−yCy on Si(001),” J. Appl. Phys., vol. 80, no. 12, pp. 6711–6715, Dec. 1996. [51] R. C. Newman and J. Wakefield, “The diffusivity of carbon in silicon,” Journal of Physics and Chemistry of Solids, vol. 19, no. 3–4, pp. 230–234, May 1961. [52] J. A. Baker, T. N. Tucker, N. E. Moyer, and R. C. Buschert, “Effect of Carbon on the Lattice Parameter of Silicon,” Journal of Applied Physics, vol. 39, no. 9, pp. 4365–4368, Aug. 1968. [53] D. Windisch and P. Becker, “Lattice distortions induced by carbon in silicon,” Philosophical Magazine A, vol. 58, no. 2, pp. 435–443, Aug. 1988. [54] G. D. Watkins and K. L. Brower, “EPR Observation of the Isolated Interstitial Carbon Atom in Silicon,” Phys. Rev. Lett., vol. 36, no. 22, pp. 1329–1332, May 1976. [55] Silicon-Germanium Carbon Alloys: Growth, Properties and Applications, 1st edition. New York: CRC Press, 2002. [56] J. Tersoff, “Enhanced Solubility of Impurities and Enhanced Diffusion near Crystal Surfaces,” Phys. Rev. Lett., vol. 74, no. 25, pp. 5080–5083, Jun. 1995. [57] R. W. Olesinski and G. J. Abbaschian, “The C−Si (Carbon-Silicon) system,” Bulletin of Alloy Phase Diagrams, vol. 5, no. 5, pp. 486–489, Oct. 1984. [58] K. Eberl, S. S. Iyer, S. Zollner, J. C. Tsang, and F. K. LeGoues, “Growth and strain compensation effects in the ternary Si1−x-yGexCy alloy system,” Applied Physics Letters, vol. 60, no. 24, pp. 3033–3035, Jun. 1992. [59] J. W. Strane, H. J. Stein, S. R. Lee, B. L. Doyle, S. T. Picraux, and J. W. Mayer, “Metastable SiGeC formation by solid phase epitaxy,” Applied Physics Letters, vol. 63, no. 20, pp. 2786–2788, Nov. 1993. [60] T. O. Mitchell, J. L. Hoyt, and J. F. Gibbons, “Substitutional carbon incorporation in epitaxial Si1−yCy layers grown by chemical vapor deposition,” Applied Physics Letters, vol. 71, no. 12, pp. 1688–1690, Sep. 1997. 109 [61] R. I. Scace and G. A. Slack, “Solubility of Carbon in Silicon and Germanium,” The Journal of Chemical Physics, vol. 30, no. 6, pp. 1551–1555, Jun. 1959. [62] J. D. Cressler, SiGe and Si Strained-Layer Epitaxy for Silicon Heterostructure Devices. CRC Press, 2010. [63] J. L. Hoyt, T. O. Mitchell, K. Rim, D. V. Singh, and J. F. Gibbons, “Comparison of Si/Si1−x-yGexCy and Si/Si1−yCy heterojunctions grown by rapid thermal chemical vapor deposition,” Thin Solid Films, vol. 321, no. 1–2, pp. 41–46, May 1998. [64] V. D. Akhmetov and V. V. Bolotov, “The effect of carbon and boron on the accumulation of vacancy-oxygen complexes in silicon,” Radiation Effects, vol. 52, no. 3–4, pp. 149–152, Jan. 1980. [65] J. P. Kalejs, L. A. Ladd, and U. Gösele, “Self‐interstitial enhanced carbon diffusion in silicon,” Applied Physics Letters, vol. 45, no. 3, pp. 268–269, Aug. 1984. [66] L. A. Ladd and J. P. Kalejs, “Self-Interstitial Injection Effects on Carbon Diffusion in Silicon At High Temperatures,” in Symposium K – Oxygen, Carbon, Hydrogen and Nitrogen in Crystalline Silicon, 1985, vol. 59. [67] U. Goesele, P. Laveant, R. Scholz, N. Engler, and P. Werner, “Diffusion Engineering by Carbon in Silicon,” MRS Online Proceedings Library, vol. 610, p. null–null, 2000. [68] R. F. Scholz, P. Werner, U. Gösele, and T. Y. Tan, “The contribution of vacancies to carbon out-diffusion in silicon,” Applied Physics Letters, vol. 74, no. 3, pp. 392–394, Jan. 1999. [69] P. Werner, U. Gösele, H.-J. Gossmann, and D. C. Jacobson, “Carbon diffusion in silicon,” Applied Physics Letters, vol. 73, no. 17, pp. 2465–2467, Oct. 1998. [70] L. D. Lanzerotti, J. C. Sturm, E. Stach, R. Hull, T. Buyuklimanli, and C. Magee, “Suppression of boron outdiffusion in SiGe HBTs by carbon incorporation,” in Electron Devices Meeting, 1996. IEDM ’96., International, 1996, pp. 249–252. [71] H. J. Osten, G. Lippert, D. Knoll, R. Barth, B. Heinemann, H. Rucker, and P. Schley, “The effect of carbon incorporation on SiGe heterobipolar transistor performance and process margin,” in Electron Devices Meeting, 1997. IEDM ’97. Technical Digest., International, 1997, pp. 803–806. [72] D. Knoll, B. Heinemann, H. J. Osten, B. Ehwald, B. Tillack, P. Schley, R. Barth, M. Matthes, K. S. Park, Y. Kim, and W. Winkler, “Si/SiGe:C heterojunction bipolar transistors in an epi-free well, single-polysilicon technology,” in Electron Devices Meeting, 1998. IEDM ’98. Technical Digest., International, 1998, pp. 703–706. [73] N. R. Zangenberg, J. Fage-Pedersen, J. L. Hansen, and A. N. Larsen, “B and P diffusion in strained and relaxed Si and SiGe,” J. Appl. Phys., vol. 94, no. 6, pp. 3883–3890, Sep. 2003. [74] J. S. Christensen, H. H. Radamson, A. Y. Kuznetsov, and B. G. Svensson, “Diffusion of phosphorus in relaxed Si1−xGex films and strained Si/Si1−xGex heterostructures,” J. Appl. Phys., vol. 94, no. 10, pp. 6533–6540, Nov. 2003. [75] H. Rucker, B. Heinemann, and A. Fox, “Half-Terahertz SiGe BiCMOS technology,” in Silicon Monolithic Integrated Circuits in RF Systems (SiRF), 2012 IEEE 12th Topical Meeting on, 2012, pp. 133 –136. [76] S. C. Jain, H. J. Osten, B. Dietrich, and H. Rucker, “Growth and properties of strained Si1−x-yGexCy layers,” Semicond. Sci. Technol., vol. 10, no. 10, p. 1289, Oct. 1995. [77] H. Rucker, B. Heinemann, W. Ropke, G. Fischer, G. Lippert, H.-J. Osten, and R. Kurps, “Modeling the effect of carbon on boron diffusion,” in , 1997 International Conference on Simulation of Semiconductor Processes and Devices, 1997. SISPAD ’97, 1997, pp. 281–284. [78] Sentaurus ProcessTM, version 2013.12, (Synopsys, Inc., Mountain View, CA, 2013), . [79] A. Sibaja-Hernandez, S. Decoutere, and H. Maes, “A comprehensive study of boron and carbon diffusion models in SiGeC heterojunction bipolar transistors,” Journal of Applied Physics, vol. 98, no. 6, p. 063530, Sep. 2005. [80] N. E. B. Cowern, P. C. Zalm, P. van der Sluis, D. J. Gravesteijn, and W. B. de Boer, “Diffusion in strained Si(Ge),” Phys. Rev. Lett., vol. 72, no. 16, pp. 2585–2588, Apr. 1994. 110 [81] T. Ito, I. Kato, T. Nozaki, T. Nakamura, and H. Ishikawa, “Plasma‐enhanced thermal nitridation of silicon,” Applied Physics Letters, vol. 38, no. 5, pp. 370–372, Mar. 1981. [82] S. T. Ahn, H. W. Kennel, J. D. Plummer, and W. A. Tiller, “Film stress‐related vacancy supersaturation in silicon under low‐pressure chemical vapor deposited silicon nitride films,” Journal of Applied Physics, vol. 64, no. 10, pp. 4914–4919, Nov. 1988. [83] M. M. Moslehi and K. C. Saraswat, “Thermal nitridation of Si and SiO2 for VLSI,” Electron Devices, IEEE Transactions on, vol. 32, no. 2, pp. 106 – 123, Feb. 1985. [84] T. Ito, T. Nozaki, H. Arakawa, and M. Shinoda, “Thermally grown silicon nitride films for high‐performance MNS devices,” Applied Physics Letters, vol. 32, no. 5, pp. 330–331, Mar. 1978. [85] S. P. Murarka, C. C. Chang, and A. C. Adams, “Thermal Nitridation of Silicon in Ammonia Gas: Composition and Oxidation Resistance of the Resulting Films,” J. Electrochem. Soc., vol. 126, no. 6, pp. 996–1003, Jun. 1979. [86] Y. Hayafuji, K. Kajiwara, and S. Usui, “Shrinkage and growth of oxidation stacking faults during thermal nitridation of silicon and oxidized silicon,” Journal of Applied Physics, vol. 53, no. 12, pp. 8639–8646, Dec. 1982. [87] J. M. Silcock and W. J. Tunstall, “Partial dislocations associated with NbC precipitation in austenitic stainless steels,” Philosophical Magazine, vol. 10, no. 105, pp. 361–389, Sep. 1964. [88] S. Mizuo, T. Kusaka, A. Shintani, M. Nanba, and H. Higuchi, “Effect of Si and SiO2 thermal nitridation on impurity diffusion and oxidation induced stacking fault size in Si,” Journal of Applied Physics, vol. 54, no. 7, pp. 3860–3866, Jul. 1983. [89] M. Karunaratne, J. Bonar, P. Ashburn, and A. Willoughby, “Suppression of boron diffusion due to carbon during rapid thermal annealing of SiGe based device materials—some comments,” Journal of Materials Science, vol. 41, no. 3, pp. 1013–1016, 2006. [90] J. M. Bonar, A. F. W. Willoughby, A. H. Dan, B. M. McGregor, W. Lerch, D. Loeffelmacher, G. A. Cooke, and M. G. Dowsett, “Antimony and boron diffusion in SiGe and Si under the influence of injected point defects,” Journal of Materials Science: Materials in Electronics, vol. 12, no. 4, pp. 219–221, 2001. [91] Y. Lin, H. Yasuda, M. Schiekofer, B. Benna, R. Wise, and G. (Maggie) Xia, “Effects of carbon on phosphorus diffusion in SiGe:C and the implications on phosphorus diffusion mechanisms,” Journal of Applied Physics, vol. 116, no. 14, p. 144904, Oct. 2014. [92] S. Kobayashi, M. Iizuka, T. Aoki, N. Mikoshiba, M. Sakuraba, T. Matsuura, and J. Murota, “Segregation and diffusion of phosphorus from doped Si1−xGex films into silicon,” Journal of Applied Physics, vol. 86, no. 10, pp. 5480–5483, Nov. 1999. [93] K. Weiser, “Theoretical calculation of distribution coefficients of impurities in germanium and silicon, heats of solid solution,” Journal of Physics and Chemistry of Solids, vol. 7, no. 2–3, pp. 118–126, Nov. 1958. [94] P. Rai-Choudhury and E. I. Salkovitz, “Doping of epitaxial silicon: Equilibrium gas phase and doping mechanism,” Journal of Crystal Growth, vol. 7, no. 3, pp. 353–360, Nov. 1970. [95] D. T. J. Hurle, R. M. Logan, and R. F. C. Farrow, “Doping of epitaxial silicon,” Journal of Crystal Growth, vol. 12, no. 1, pp. 73–75, Jan. 1972. [96] W. H. Shepherd, “Doping of Epitaxial Silicon Films,” J. Electrochem. Soc., vol. 115, no. 5, pp. 541–545, May 1968. [97] C. D. Thurmond and J. D. Struthers, “Equalibrium Thermochemistry of Solid and Liquid Alloys of Germanium and of Silicon. II. The Retrograde Solid Solubilities of Sb in Ge, Cu in Ge and Cu in Si,” J. Phys. Chem., vol. 57, no. 8, pp. 831–835, Aug. 1953. [98] S. M. Hu, D. C. Ahlgren, P. A. Ronsheim, and J. O. Chu, “Experimental study of diffusion and segregation in a Si-(GexSi1−x) heterostructure,” Phys. Rev. Lett., vol. 67, no. 11, pp. 1450–1453, Sep. 1991. 111 [99] S. Kobayashi, T. Aoki, N. Mikoshiba, M. Sakuraba, T. Matsuura, and J. Murota, “Segregation and diffusion of impurities from doped Si1−xGex films into silicon,” Thin Solid Films, vol. 369, no. 1–2, pp. 222–225, Jul. 2000. [100] R. F. Lever, J. M. Bonar, and A. F. W. Willoughby, “Boron diffusion across silicon–silicon germanium boundaries,” Journal of Applied Physics, vol. 83, no. 4, pp. 1988–1994, Feb. 1998. [101] S. M. Hu, “General Theory of Impurity Diffusion in Semiconductors via the Vacancy Mechanism,” Phys. Rev., vol. 180, no. 3, pp. 773–784, Apr. 1969. [102] S. M. Hu, “Diffusion and segregation in inhomogeneous media and the Ge_{x}Si_{1-x} heterostructure,” Phys. Rev. Lett., vol. 63, no. 22, pp. 2492–2495, Nov. 1989. [103] C. Ahn and S. T. Dunham, “Calculation of dopant segregation ratios at semiconductor interfaces,” Phys. Rev. B, vol. 78, no. 19, p. 195303, Nov. 2008. [104] H.-M. You, U. M. Gösele, and T. Y. Tan, “Simulation of the transient indiffusion‐segregation process of triply negatively charged Ga vacancies in GaAs and AlAs/GaAs superlattices,” Journal of Applied Physics, vol. 74, no. 4, p. 2461, Aug. 1993. [105] T. Y. Tan, “Mass transport equations unifying descriptions of isothermal diffusion, thermomigration, segregation, and position-dependent diffusivity,” Applied Physics Letters, vol. 73, no. 18, p. 2678, Nov. 1998. [106] T. Y. Tan, R. Gafiteanu, H.-M. You, and U. Gösele, “Diffusion‐segregation equation for modeling in heterostructures,” AIP Conference Proceedings, vol. 306, no. 1, p. 478, Jun. 1994. [107] Paul G. Shewmon, “Problem 4-8 in Chapter 4,” in Diffusion in Solids, p. 136. [108] M. E. Glicksman, Diffusion in Solids: Field Theory, Solid-State Principles, and Applications. Wiley, 1999. [109] P. W. Atkins, Physical Chemistry. W.H. Freeman, 1994. [110] L. J. Giling and J. Bloem, “The incorporation of phosphorus in silicon; The temperature dependence of the segregation coefficient,” Journal of Crystal Growth, vol. 31, pp. 317–322, Dec. 1975. [111] MediciTM, version 2010.03, Synopsys, Inc., Mountain View, CA, 2010, . [112] R. Kube, H. Bracht, E. Hüger, H. Schmidt, J. L. Hansen, A. N. Larsen, J. W. Ager III, E. E. Haller, T. Geue, and J. Stahn, “Contributions of vacancies and self-interstitials to self-diffusion in silicon under thermal equilibrium and nonequilibrium conditions,” Phys. Rev. B, vol. 88, no. 8, p. 085206, Aug. 2013. 112 Appendix A Numerical fitting with TCAD softwares The simulation work for diffusivity extraction in chapter 4 and chapter 5 was performed using TSUPREM-4 software. The established carbon/Ge model was implemented both in TSUPREM-4 and Sprocess, which is another multidimensional process modeling environment. Example codes for each of the above purpose were provided as below. Appendix A1 P diffusivity ratio extraction in TSUPREM-4 $ This program takes as-grown SIMS data as input, simulates and compares annealed data vs. model prediction $ Set the filenames and annealing parameters assign name=T n.val=1000 assign name=AnnealTime n.val=30/60 assign name=Gefactor n.val=1.5 assign name=Pfactor n.val=0.6 assign name=GeInput c.val=KA-R_Ge.txt assign name=GeAnnealed c.val=KCA_Ge.txt assign name=PInput c.val=KA-R_P.txt assign name=PAnnealed c.val=KCA_P.txt assign name=XgeSubstrate n.val=0.00 $ Define dopant diffusivity ratio IMPURITY IMP=PHOSPH MAT=SILICON MODEL=PinSiGe DI.FAC=@Pfactor METHOD MODEL=PinSige ENABLE 113 $ Define mesh line x loc=0 spacing=0.001 line x loc=0.01 spacing=0.01 line y loc=0 spacing=0.001 line y loc=0.06 spacing=0.0002 line y loc=0.16 spacing=0.0002 line y loc=0.2 spacing=0.0010 init boron=1e15 material=silicon method st.histo ^skip.sil viscoela method pd.full method ^dif.adap $------------------------------------ select z=1 title='mesh' plot.2d grid x.min=-0.01 x.max=0.01 y.min=0 y.max=0.20 c.grid=2 $ ---- load and plot initial profile profile impurity=Ge inf=@GeInput select z=Ge title='DrDc modeling T=@T C @AnnealTime mins ' plot.1d x.min=0.05 x.max=0.18 y.min=0 y.max=1e22 $label x=0.015 y=0.6e22 label="b1=8.1 b2=23 Xge0=@XgeSubstrate" 114 PROFILE PHOSPH IN.FILE=@PInput SELECT Z=LOG10(PHOSPHORUS) PLOT.1D ^CLEAR x.min=0.05 x.max=0.18 y.min=16 y.max=20 LINE.TYP=1 COLOR=1 $ ---- here is the diffusion process step $diffusion temp=400 t.final=@T time=(@T-400)/100/60 inert diffusion temp=@T time=@AnnealTime inert $diffusion temp=@T t.final=400 time=(@T-400)/100/60 inert $ ---- plot the diffused profile on top of the initial profile SELECT Z=Ge PLOT.1D ^CLEAR x.min=0.05 x.max=0.18 y.min=0 y.max=1e22 LINE.TYP=1 COLOR=2 SELECT Z=LOG10(PHOSPHORUS) PLOT.1D ^CLEAR x.min=0.05 x.max=0.18 y.min=16 y.max=20 LINE.TYP=1 COLOR=2 label x=0.13 y=19.5 label="DrDc model prediction" label x=0.13 y=19.0 label="Ge Diffusivity scale: X@Gefactor" label x=0.13 y=18.8 label="P Diffusivity scale: X@Pfactor" label x=0.13 y=18.5 label="T = @T C" label x=0.13 y=18.3 label="Time = @AnnealTime mins" 115 $ ---- output the simulation results in a txt file SELECT Z=Ge EXTRACT OUT.FILE=annealed_simulation.txt PREFIX="% Ge profile simulated" FOREACH DEPTH (0 TO 0.08 STEP 0.0002) EXTRACT SILICON X=0.0 DISTANCE=@{DEPTH} Y.EXT VAL.EXT END EXTRACT CLOSE SELECT Z=Ge/5e20 PRINT.1D X.VAL=0 X.MAX=0.2 OUT.FILE=KA-R-@T-@AnnealTime-min-Gefactor-@Gefactor-pd.full.txt SELECT Z=PHOSPHORUS PRINT.1D X.VAL=0 X.MAX=0.2 OUT.FILE=KA-R-@T-@AnnealTime-min-Pfactor-@Pfactor-pd.full.txt impurity name=verticalline new assign name=reference c.val=vertical-line.txt PROFILE IMPURITY=verticalline IN.FILE=@reference replace SELECT Z=LOG10(verticalline) PLOT.1D ^CLEAR x.min=0.05 x.max=0.18 y.min=16 y.max=20 LINE.TYP=1 COLOR=2 $ ---- load SIMS data of annealed sample, and overlay with initial and simulation profiles PROFILE IMPURITY=Ge IN.FILE=@GeAnnealed replace 116 SELECT Z=Ge PLOT.1D ^CLEAR x.min=0.05 x.max=0.18 y.min=0 y.max=1e22 LINE.TYP=1 COLOR=4 PROFILE PHOSPHORUS IN.FILE=@PAnnealed replace SELECT Z=LOG10(PHOSPHORUS) PLOT.1D ^CLEAR x.min=0.05 x.max=0.18 y.min=16 y.max=20 LINE.TYP=1 COLOR=4 label x=0.14 y=17.5 label="Annealed SIMS" c.line=4 line.typ=1 label x=0.013 y=.4e22 label="SIMS data annealed" c.line=4 line.typ=1 Appendix A2 Carbon/Ge model implementation in TSUPREM-4 $ Set the filenames and annealing parameters assign name=T n.val=1000 assign name=AnnealTime n.val=30/60 assign name=Gefactor n.val=1.5 assign name=Cfactor n.val=0.08 assign name=XgeSubstrate n.val=0.00 assign name=GeInput c.val=K_asgrown_Ge.txt assign name=PInput c.val=K_asgrown_P.txt assign name=CInput c.val=K_asgrown_C.txt assign name=GeAnnealed c.val=KCA_30s_Ge.txt assign name=PAnnealed c.val=KCA_30s_P.txt assign name=CAnnealed c.val=KCA_30s_C.txt 117 $-------------define mesh ----------------------- line x loc=0 spacing=0.001 line x loc=0.01 spacing=0.01 line y loc=0 spacing=0.001 line y loc=0.04 spacing=0.0002 line y loc=0.20 spacing=0.0002 line y loc=0.23 spacing=0.001 init boron=1e15 material=silicon method st.histo ^skip.sil viscoela method pd.full method ^dif.adap $-------------------------------------------------- select z=1 title='mesh' plot.2d grid x.min=0 x.max=0.01 y.min=0 y.max=0.20 c.grid=2 $--------------------------------------- $ Define Ge as new impurity and method impurity name=Ge new method variable=Ge none abs.err=0.1 118 profile impurity=Ge inf=@GeInput PLOT.1D x.min=0.00 x.max=0.25 y.min=0 y.max=1.5e22 LINE.TYP=1 COLOR=1 $ ---- define interdiffusion model and parameters intermed name=Xge express=Ge/5e22 intermed name=Xge0 value=@XgeSubstrate intermed name=b1 value=8.1 intermed name=b2 value=23 intermed name=D00 value=310*@Gefactor intermed name=Eav value=4.66 intermed name=D01 express=D00*exp(-Eav/kT) intermed name=D02 value=D01*exp(b1*Xge0) intermed name=D express=(Xge<=Xge0 ? D01*exp(b1*Xge):D02*exp(b2*(Xge-Xge0))) equation variable=Ge mat=Si addtoexp=DIV(D*GRAD(Ge)) $------------------------------------------ $ Define C Model $impurity name=Cnew new $method variable=Cnew none abs.err=0.1 PROFILE impurity=carbon IN.FILE=@CInput SELECT Z=LOG10(carbon) 119 PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=16 y.max=20 LINE.TYP=1 COLOR=1 IMPURITY IMP=CARBON MAT=SILICON MODEL=Cscale DI.FAC=@Cfactor METHOD MODEL=Cscale ENABLE $ Define P Model $ FIP = 1-0.02*Xge+xGe^2-13*xGe^3 intermed name=XCarbon express=carbon/5e22 intermed name=FIP express=1-0.02*Xge+xGe^2-13*xGe^3 intermed name=FVP express=1-FIP intermed name=Pfactor express=FIP/(1+XCarbon*5500)+FVP*(1+XCarbon*5500) IMPURITY IMP=PHOSPH MAT=SILICON MODEL=PinSiGeC DI.FAC=Pfactor METHOD MODEL=PinSiGeC ENABLE $ ---- load and plot initial profile PROFILE PHOSPH IN.FILE=@PInput SELECT Z=LOG10(PHOSPHORUS) PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=16 y.max=20 LINE.TYP=1 COLOR=1 $ ---- here is the diffusion process step $diffusion temp=400 t.final=@T time=(@T-400)/100/60 inert diffusion temp=@T time=@AnnealTime inert $diffusion temp=@T t.final=400 time=(@T-400)/100/60 inert 120 $ ---- plot the diffused profile on top of the initial profile SELECT Z=Ge PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=0 y.max=1.5e22 LINE.TYP=1 COLOR=2 SELECT Z=LOG10(PHOSPHORUS) PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=16 y.max=20 LINE.TYP=1 COLOR=2 SELECT Z=LOG10(CARBON) PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=16 y.max=20 LINE.TYP=1 COLOR=2 label x=0.19 y=19.7 label="T = @T C" label x=0.19 y=19.4 label="Time = @AnnealTime mins" label x=0.19 y=19.0 label="DrDc X @Gefactor" label x=0.19 y=18.7 label="Cscale = @Cfactor" SELECT Z=PHOSPHORUS PRINT.1D X.VAL=0 X.MAX=0.22 OUT.FILE=UBC_KCA_P_@T-@AnnealTime-min-CDIFAC=@Cfactor-Gefactor=@Gefactor-pd.full.txt SELECT Z=CARBON PRINT.1D X.VAL=0 X.MAX=0.22 OUT.FILE=UBC_KCA_C_@T-@AnnealTime-min-CDIFAC=@Cfactor-Gefactor=@Gefactor-pd.full.txt SELECT Z=Ge/5e20 121 PRINT.1D X.VAL=0 X.MAX=0.22 OUT.FILE=UBC_KCA_Ge_@T-@AnnealTime-min-CDIFAC=@Cfactor-Gefactor=@Gefactor-pd.full.txt $ ---- load SIMS data of annealed sample, and overlay with initial and simulation profiles PROFILE IMPURITY=Ge IN.FILE=@GeAnnealed replace SELECT Z=Ge PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=0 y.max=1.5e22 LINE.TYP=1 COLOR=4 PROFILE PHOSPHORUS IN.FILE=@PAnnealed replace SELECT Z=LOG10(PHOSPHORUS) PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=16 y.max=20 LINE.TYP=1 COLOR=4 PROFILE IMPURITY=CARBON IN.FILE=@CAnnealed replace SELECT Z=LOG10(CARBON) PLOT.1D ^CLEAR x.min=0.00 x.max=0.25 y.min=16 y.max=20 LINE.TYP=1 COLOR=4 Appendix A3 Carbon/Ge model implementation in Sentaurus Process #header set DEBUG 1 pdbSet ImplantData DoseControl WaferDose #endheader pdbSet Silicon Dopant DiffModel pair 122 math coord.ucs #----------------------------------------------------------------------# mgoals on min.normal.size= 5 accuracy= 1 \ normal.growth.ratio= 1.5 refinebox name= Interface \ min.normal.size= 1 \ interface.materials= "Silicon" \ min= "-0.01 0.2" max= "+0.01 1.0" refinebox name= All min= "-1 0" max= "1 1" \ xrefine= "0.01 0.005 0.01" yrefine= "0.1" add all #----------------------------------------------------------------------# line x location= -0.2 spacing= 0.05 line x location= 0 spacing= 0.005 tag= sitop line x location= 0.1 spacing= 0.05 tag= sibot line y location= 0 spacing= 0.5 tag= left line y location= 1 spacing= 0.5 tag= right region substrate silicon xlo= sitop xhi= sibot ylo= left yhi= right init field= Boron concentration= 1e15 !DelayFullD #UserAddEqnTerm Silicon Boron "{1*Boron/2*10*grad(Germanium)}" #UserAddEqnTerm Silicon Boron "{1e-13*grad(Boron)}" 123 #Term name= DiffFactor add Silicon eqn = "exp(2)" term name=BoronIntDiffFactor add Silicon eqn = "5*exp(2)" term name=BoronVacDiffFactor add Silicon eqn = "5" #pdbset silicon boron BoronDiffFactor #MultiplyTerm Si BoronDiffFactor eqn="5" #MultiplyTerm Silicon BoronInt "exp(2)" #----------------------------------------------------------------------# #split @SiGe_C@ refinebox name= SiGe min= "-0.1 0" max= "0 1" \ xrefine= "0.001" yrefine= "0.05" add all strain_profile species= Germanium silicon ratio= {0 1} strain= {0 0.0425} pdbSetBoolean Silicon Mechanics UpdateStrain 1 deposit material= {Silicon} type= isotropic rate= {1.0} time= 0.2 #-Setting Boron profile in newly deposited SiGe layer sel z= "-0.055 < x && x < -0.045 ? 1e20 : 1e12" name= NewP Silicon sel z= "Boron + NewP" name= Boron store sel name= NewP delete #-Setting Germanium profile in newly deposited SiGe layer sel z= "-0.08 < x && x < -0.02 ? 1.0e22 : 1e10" \ 124 name= Germanium Silicon store if { $DEBUG } { SetPlxList { Boron Germanium } WritePlx n@node@_SiGe_as_grown.plx y= 0.5 } #----------------------------------------------------------------------# #split @RTA@ diffuse time=0.5 temp=1000 if { $DEBUG } { SetPlxList {PTotal Boron BoronInt BoronVac int interstitial vac vacancy Germanium} WritePlx n@node@_SiGe_RTA.plx y=0.5 } struct tdr= n@node@ 125 Appendix B Numerical Fitting with MATLAB Appendix B1 Numerical solver for diffusion-segregation equation In chapter 6, we re-derived a coupled model (Eq. 6.18) to simulate the diffusion-segregation phenomenon which occurred in inhomogeneous material system. To numerically simulate with this model, we developed a set of MATLAB codes. Below is an example: close all; clear all; clc; Temp = 900; % Unit: celcius degree. Anneal temperature. k_eV = 8.617E-5 ; % eV/K kT = k_eV*(Temp+273.15); time = 1200.0; % Unit: second. Anneal time. xrange = 0.15; % Unit: um dx = 0.0002; % Unit: um. dt = 0.02; % Unit: second. Time step. numx = int16(1+xrange/dx); % number of grid points in x=profiledepth/step numt = time/dt; % number of time steps to be iterated over x = 0:dx:xrange; % vector of x values, to be used for plotting P_input=importdata('S5A_P.txt'); P_interpolated=SIMS_Interpolation(P_input,xrange,dx); Ge_input=importdata('S5C_Ge.txt'); Ge_interpolated=SIMS_Interpolation(Ge_input,xrange,dx); % Input profile. Unit is micron. Step is 0.0001 um. Ge = zeros(numx,numt); %initialize Ge to zero P = zeros(numx,numt); %initialize P to zero Pdose = zeros(1,numt); % define diffusion and segregation parameters 126 %DP = 2.9e-14; %DP = 1.5E-15; DrDcFactor = 0.3; Dscale = 1.8; EGe = 0.6; t(1) = 0; % at t=0 Ge(:,1) = Ge_interpolated(:,2); P(:,1) = P_interpolated(:,2); %iterate difference equation - note that Ge(1,j) and Ge(numx,j) always remain 0 for j=1:numt-1 t(j+1) = t(j) + dt; for i=2:numx-1 DGe = DrDcFactor*310*exp(-4.66/8.617e-5/(273.15+Temp))*exp(23*Ge(i,j)/5e22); GeTerm = (dt/dx^2)*(Ge(i+1,j) - 2*Ge(i,j) + Ge(i-1,j)); Ge(i,j+1) = Ge(i,j) + DGe*1e8*GeTerm; %PTerm = (dt/dx^2)*(P(i+1,j) - 2*P(i,j) + P(i-1,j)); PRight = (P(i+1,j)+P(i,j))/2; PLeft = (P(i-1,j)+P(i,j))/2; xGeR = (Ge(i+1,j)+Ge(i,j))/2/5e22; xGeL = (Ge(i-1,j)+Ge(i,j))/2/5e22; fGe = EGe/(kT*5E22); ni_Si = 3.87e16*exp(-0.605/kT+1.5*log(Temp+273)); % 3.87e16*exp(-0.605/kt+1.5*log(Temp)) 127 dEgR = 0.835*xGeR^2-1.01*xGeR; dEgL = 0.835*xGeL^2-1.01*xGeL; ni_SiGeR= ni_Si*exp(-dEgR/(2*kT)); % Ni*exp(GeFrac*dEg/(2*kt)) from TS4 ni_SiGeL= ni_Si*exp(-dEgL/(2*kT)); DIX = 3.85*exp(-3.66/kT); DIM = 4.44*exp(-4/kT); DIMM = 44.2*exp(-4.37/kT); DPfermiR = DIX+DIM*PRight/ni_SiGeR+DIMM*(PRight/ni_SiGeR)^2; DPfermiL = DIX+DIM*PLeft/ni_SiGeL+DIMM*(PLeft/ni_SiGeL)^2; DPR = Dscale*DPfermiR*(1+xGeR); DPL = Dscale*DPfermiL*(1+xGeL); DiffR = DPR*(dt/dx^2)*(P(i+1,j) - P(i,j)); DiffL = DPL*(dt/dx^2)*(P(i-1,j) - P(i,j)); SegR = DPR*(dt/dx^2)*PRight*fGe*(Ge(i+1,j) - Ge(i,j)); SegL = DPL*(dt/dx^2)*PLeft*fGe*(Ge(i-1,j) - Ge(i,j)); PsegGe = SegR + SegL; P(i,j+1) = P(i,j) + 1e8*(DiffR + DiffL) + 1e8*PsegGe; end P(1, j+1) = P(2,j+1); P(numx,j+1) = P(numx-1,j+1); end P_Seg_final = P(:,numt); P_input=importdata('S5B_P.txt'); P_annealed=SIMS_Interpolation(P_input,xrange,dx); P_Seg_S5B_EGe = P(:,numt); Figure; % plot profiles %semilogy(x,Ge(:,1),'color','black'); %hold on; %semilogy(x,Ge(:,numt/10),'color','green'); %semilogy(x,Ge(:,numt/2),'color','blue'); %semilogy(x,Ge(:,numt),'Color','red'); 128 semilogy(x,P(:,1),'color','black', 'LineWidth', 2); hold on; semilogy(x,P(:,numt*1/10),'color','red', 'LineWidth', 2); semilogy(x,P(:,numt*2/10),'color','green', 'LineWidth', 2); semilogy(x,P(:,numt*3/10),'color','blue', 'LineWidth', 2); semilogy(x,P(:,numt*4/10),'color','cyan', 'LineWidth', 2); semilogy(x,P(:,numt*5/10),'color','magenta', 'LineWidth', 2); semilogy(x,P(:,numt*6/10),'color','yellow', 'LineWidth', 2); semilogy(x,P(:,numt*7/10),'color','red', 'LineWidth', 2); semilogy(x,P(:,numt*8/10),'color','green', 'LineWidth', 2); semilogy(x,P(:,numt*9/10),'color','blue', 'LineWidth', 2); semilogy(x,P(:,numt),'Color','black', 'LineWidth', 2); P_anneal=importdata('S5B_P.txt'); %Ge_anneal=importdata('S4B_Ge.txt'); %semilogy(Ge_anneal(:,1),Ge_anneal(:,2)); semilogy(P_anneal(:,1),P_anneal(:,2)); xlabel('Depth (um)'); ylabel('Concentration (atom/cm3)'); % monitor dose variation Pdose = sum(P)*dx*1e-4; Figure; plot(t,Pdose); xlabel('time (second)'); ylabel('P dose (cm-2)'); Appendix B2 Numerical search for best fitting with least-square method The fitting quality of a simulation is represented by how well the simulated profiles match to the experimental data. For this purpose, we developed a code to find the best fitting condition with least-square method. Below is an example of the code: 129 close all; clear all; clc; % S = sigma(ri^2) = (A1-f1)^2 +(As-f2)^2+(A3-f3)^2 % fi = FI / (1+k*Ci)+(1-FI)*(1+k*Ci) A = [0.8;0.45;0.5]; B = [0.55;0.53;0.75]; C = [0.058;0.093;0.18]; plot(C,A,C,B); dx = 0.01; FI_limit = 0.75; FI = FI_limit:dx:1; %FI = 0.91; K = 1:1:50; Num_FI = length(FI); Num_K = length(K); S = zeros(Num_K,Num_FI); for i=1:Num_K for j=1:Num_FI f1 = FI(1,j)/(1+K(1,i)*C(1,1))+(1-FI(1,j))*(1+K(1,i)*C(1,1)); rA1 = A(1,1)-f1; rB1 = B(1,1)-f1; f2 = FI(1,j)/(1+K(1,i)*C(2,1))+(1-FI(1,j))*(1+K(1,i)*C(2,1)); rA2 = A(2,1)-f2; rB2 = B(2,1)-f2; f3 = FI(1,j)/(1+K(1,i)*C(3,1))+(1-FI(1,j))*(1+K(1,i)*C(3,1)); rA3 = A(3,1)-f3; rB3 = B(3,1)-f3; S(i,j) = rA1^2+rB1^2+rA2^2+rB2^2+rA3^2+rB3^2; end end 130 Appendix C Remarks from the author Let’s recall a question: what is a material? In the engineering domain, it usually refers to physical objects which own certain utilities. In this sense, it should refer to physical reality and existence in our world (i.e., here we exclude virtual forms of material such as “big data”). Natural resources become materials as soon as they are harnessed by another part of nature: the living organisms, especially the human species. The human history of the utilization, renovation and invention of materials can be dated back beyond written records. If we parallelize the chronological evolvement of materials and the human understanding of them, one can immediately realize their strong interaction and mutual impacts. In the Stone Age, one of the best things human beings could do is to “shape” the rocks found from nature, which may be the earliest form of materials engineering. We may describe it as a “three-dimensional geometrical art”. Then what is the fourth dimension in materials engineering? Perhaps it should be attributed to heat, or thermal processing, which prevailed in the Bronze and the Iron Age. Here, we should note that the concept of heat as perceived by us today is significantly different from the people at that time when these technologies were invented. Since the ancient people found that heat can be transferred from one body to another, they hypothesized that heat may be a form of substance. Today, we interpret it as a form of energy. Is this difference important? The answer is yes, since it is usually our understanding and beliefs which lead us towards daily practices. Then, how did we arrive this refreshed understanding of heat and many things beyond? For example, the atomic theory was established about one century ago and the general theory of relativity was published 100 years ago, and they fundamentally changed our view of nature and universe. Perhaps we all need to thank Sir Isaac Newton, who established a “methodology” for the quest of nature which he labeled as “experimental philosophy”, or “natural philosophy” as revised by himself in later edition of his “Mathematical Principle of Natural Philosophy”. Today, we call it “science”. Based on the works of many people before him, he laid out the principles of experimentation which guided the later generations in the investigation of nature, the universe and human himself. Since then, tremendous progress has been achieved in regards to the expansion of experimentally-tested knowledge and its impact on engineering, technology, medicine, and even social aspects such as law and religion, to name just a few. Now, let’s finish this thesis with two quotes: 131 — A.H. Compton, Journal of Applied Physics 1, 13 (1931). — A Cree proverb engraved on the ground outside of UBC Law School