UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 1 ENHANCEMENTS IN CMOS DEVICE SIMULATION FOR SINGLE EVENT EFFECTS By DANIEL J. CUMMINGS A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2010 PAGE 2 2 2010 Daniel J. Cummings PAGE 3 3 To Julie PAGE 4 4 ACKNOWLEDGMENTS First and foremost, I would like to thank Dr. Mark E. Law This work would not have been possible without his support and guidance. My work experience in graduate school was incredibly enjoyable and this was due to having him as my advisor. Next, I would like to thank Hyunwoo Park for his invaluable insight, patience and consistent hard work on many of the projects that we shared. Additio nally, I want to thank Srivatsan for all of the useful discussions on strained silicon technology and MOSFET physics. I want to acknowledge Saurabh and Nicole for their assistance with FLOODS code and appreciate all the other SWAMP he privilege of working with along the way I am very grateful to Dr. Scott E. Thompson, Dr. Jing Guo and Dr. Brent Gila for their support and guidance on my supervisory committee. I want to thank Steve Cea and Tom Linton from Intel for their support on t he discretization work. From Vanderbilt, I would like to thank Dr. Schrimpf and Dr. Witulski for their input on the bulk mobility modeling. Also, I want to acknowledge the people at the Naval research lab such as Charles Beckler, Philomena West, Dan Crute and many others for the invaluable experience I gained while working there during my undergraduate years. I would especially like to thank my parents Lois and Joseph, for their love and encouragement I want to acknowledge my siblings Stacie, Jonathan and Andrew for their support and friendship. My group of friends in Gainesville was always there to balance all of the hard schoolwork with fun. I want to thank Timothy, Jason, John, Rich, Ricky and so many others for some great memories. Last, but certainly not least, this work would not have been possible without my wife Julie. I want to thank her for her overwhelming love and support in my efforts. PAGE 5 5 TABLE OF CONTENTS page ACKNOWLEDGMENTS ................................ ................................ ................................ .. 4 LIST OF TABLES ................................ ................................ ................................ ............ 9 LIST OF FIGURES ................................ ................................ ................................ ........ 10 LIST OF ABBREVIATIONS ................................ ................................ ........................... 18 ABSTRACT ................................ ................................ ................................ ................... 19 CHAPTER 1 INTRODUCTION AND BACKGROUND ................................ ................................ 22 1.1 Motivation ................................ ................................ ................................ ...... 22 1.2 Brief Overview of Single event Effects ................................ .......................... 23 1.2.1 Brief History of Single event Effects ................................ .................... 24 1.2.2 Radiation Sources ................................ ................................ ............... 25 1.2.3 Example: Single event Upset in a 6T SRAM ................................ ....... 26 1.3 CMOS Scaling and Susceptibility ................................ ................................ .. 27 1.4 Single Event Device Simulation Challenges ................................ ................. 29 1.5 FLOODS Simulation Tool ................................ ................................ .............. 30 2 PHYSICAL MECHANISMS OF SINGLE EVENT EFFECTS ................................ .. 37 2.1 Introduction ................................ ................................ ................................ ... 37 2.2 Carrier Generation ................................ ................................ ........................ 37 2.3 Particle Strike Models ................................ ................................ ................... 39 2.3.1 Linear Energy Transfer ................................ ................................ ....... 40 2.3.2 Heavy Ion Modeling ................................ ................................ ............ 41 2.3.3 Pulsed Laser Modeling ................................ ................................ ....... 42 2.4 Charge Collection Mechanisms ................................ ................................ ..... 44 2 .4.1 Baseline Simulation Structure ................................ ............................. 44 2.4.2 The Basics of Charge Transport ................................ ......................... 45 2.4.3 Analytic Approximations ................................ ................................ ...... 46 2.4.4 Doping Profile Effects ................................ ................................ .......... 47 2.4.5 Mobility ................................ ................................ ................................ 49 2.4.6 Charge Conservation ................................ ................................ .......... 50 2.4.7 Recombination ................................ ................................ .................... 52 2.4.7.1 Auger Recombination ................................ ................................ .. 53 2.4.7.2 SRH Recombination ................................ ................................ .... 53 2.4.7.3 The Impact of Recombination on Charge Collection .................... 54 2.4.8 Bandgap Narrowing ................................ ................................ ............ 55 2.5 Summary ................................ ................................ ................................ ....... 56 PAGE 6 6 3 DISCRETIZATION METHODS FOR SEE SIMULATIONS ................................ ..... 75 3.1 Introduction ................................ ................................ ................................ ... 75 3.2 Discretization Overview ................................ ................................ ................. 76 3.2.1 Finite Volume Discretization ................................ ............................... 78 3.2.2 Finite Element Discretization ................................ .............................. 79 3.3 Simulation Methodology ................................ ................................ ................ 80 3.4 Simulation Results ................................ ................................ ........................ 81 3.4.1 Short Channel MOSFET results ................................ .......................... 81 3.4.2 Charge Collection simulations ................................ ............................. 82 3.5 Discretization M ethod Summary ................................ ................................ ... 84 4 DEVICE GRID AND BOUNDARY SCHEMES ................................ ........................ 91 4.1 Introduction ................................ ................................ ................................ ... 91 4.2 Adaptive Gridding ................................ ................................ .......................... 91 4.2.1 Minimizing Discretization Error ................................ ............................ 92 4.2.2 Simulation Time Tradeoff ................................ ................................ .... 93 4.2.3 Adaptive Grid Scheme Methodology ................................ ................... 94 4.2.4 Simulation Results ................................ ................................ .............. 97 4.2.4.4 N+/ P Diode Simulation ................................ ................................ 98 4.2.4.5 NMOS Simulation ................................ ................................ ...... 100 4.2.5 Adaptive Grid Summary ................................ ................................ .... 101 4.3 Boundary Sinks ................................ ................................ ........................... 102 4.3.1 Boundary Condition Overview ................................ ........................... 102 4.3.2 Simulation Results ................................ ................................ ............ 105 4.3.3 Boundary Sink Summary ................................ ................................ .. 107 5 IMPACT OF STRAINED SILICON ON SINGLE EVENT EFFECTS ..................... 118 5.1 Motivation and Background ................................ ................................ ......... 118 5.2 Brief Overview of the Physics of Strained Silicon ................................ ........ 119 5.3 Piezoresistance Mobility Mo del ................................ ................................ ... 121 5.3.1 Linear Elasticity ................................ ................................ ................. 121 5.3.2 The Strain and Stress Tensors ................................ ......................... 123 5.3.3 Piezoresistance Definition ................................ ................................ 126 5.3.4 Piezoresistive Mobility Model Implementation ................................ ... 130 5.4 Uniaxial Strained Si Diode ................................ ................................ .......... 133 5.4.1 Experimental Setup ................................ ................................ ........... 134 5.4.2 Comparison of Experimental and Simulation Results ....................... 135 5.4.3 Uniaxially Strained Si Diode Summary ................................ ............. 136 5.5 Predictions for Strained Si MOSFETs ................................ ......................... 137 5 .5.1 Simulation Setup Overview ................................ ............................... 138 5.5.2 NMOS Simulation Results ................................ ................................ 140 5.5.3 PMOS Simulation Results ................................ ................................ 142 5.5.4 Impact of STI on Single event Transients ................................ ......... 143 5.5.5 STI Simulation Results ................................ ................................ ...... 144 PAGE 7 7 5.5.6 Strained Si MOSFET Summary ................................ ........................ 146 6 BULK MOBILITY MODELING FOR SINGLE EVENT EFFECTS .......................... 174 6.1 Introduction ................................ ................................ ................................ 174 6.2 Overview of Existing Bulk Mobility Models ................................ .................. 176 6.3 High Injection Mobility Model ................................ ................................ ...... 178 6.3.1 Majority Carrier Modeling ................................ ................................ .. 179 6.3.2 Minority Carrier Modeling ................................ ................................ .. 181 6.3.3 Carrier Carrier Scattering ................................ ................................ .. 183 6.3.4 Simulation Results and Discussion ................................ ................... 184 6.3.5 Experiment Setup ................................ ................................ .............. 185 6.3.6 Gen erated Carrier Distribution ................................ .......................... 186 6.3.7 Simulation Set 1 Results Experimental Comparison ...................... 186 6.3.8 Simulation Set 2 Results Io n Strike ................................ ................ 187 6.3.9 Simulation Set 3 Results Epitaxial Diode ................................ ....... 187 6.3.10 Computational Comparison ................................ ............................... 188 6.3.11 Summary ................................ ................................ ........................... 188 6.4 General Purpose Mobility Model ................................ ................................ 189 6.4.1 Lattice Scattering ................................ ................................ .............. 191 6.4.2 Majority Impurity Scattering ................................ ............................... 191 6.4.3 Minority Impurity Scattering and Charge Screening .......................... 193 6.4.4 Electron Hole Scattering and Charge Screening .............................. 194 6.4.5 Temperature Dependence ................................ ................................ 197 6.4.6 S imulation Results ................................ ................................ ............ 197 6.4.6.6 Bipolar N/P/N transistor simulation ................................ ............ 198 6.4.6.7 N+/P diode simulation ................................ ................................ 200 6.4.7 Computational Comparison ................................ ............................... 201 6.4.8 Summary ................................ ................................ ........................... 201 6.5 Interface Mobility Models ................................ ................................ ............. 202 6.5.1 Lombardi Model ................................ ................................ ................ 202 6.5.2 Velocity Saturation Model ................................ ................................ 203 7 SU MMARY, CONCLUSIONS AND RECOMMENTATIONS FOR FUTURE WORK ................................ ................................ ................................ ................... 225 7.1 Summary and Conclusions ................................ ................................ ......... 225 7.2 Recommendations for Future Work Equation Chapter (Next) Section 1 ...... 227 7.2.1 Carrier Ge neration with Hydrodynamic Transport ............................. 228 7.2.2 Bandgap Narrowing ................................ ................................ .......... 229 7.2.3 3 D Adaptive Gridding ................................ ................................ ....... 229 7.2.4 Single Event Experiments ................................ ................................ 230 APPENDIX DERIVATION OF TRANSFORMABLE PIEZOCOEFFICIENTS ............. 231 LIST OF REFERENCES ................................ ................................ ............................. 236 PAGE 8 8 BIOGRAPHICAL SKETCH ................................ ................................ .......................... 245 PAGE 9 9 LIST OF TABLES Table page 1 1 Sing le Event Effects terminology description. ................................ .................... 31 1 2 Spacecraft for which single event effects have impacted ................................ ... 31 2 1 Parameters for c alculating LET for various semiconductor materials. ................ 57 2 2 Experimental parameters for the single photon absorption pulsed laser ............ 57 2 3 Standard coefficients for Auger recombination model ................................ ........ 57 2 4 Standard coefficien ts for SRH recombination model ................................ .......... 58 2 5 Para meters for silicon bandgap narrowing models. ................................ ............ 5 8 3 1 Physical model sets used for the simulation comparisons ................................ .. 85 4 1 Overview of t he adaptive grid simulation variables ................................ ........... 108 5 1 Values of piezoresistance ( ) coefficients ................................ ........................ 147 6 1 Majority Carrier Mobility Fitti ng Parameters at 300 K. ................................ ...... 204 6 2 Temperature Dependence Fitting Parameters ................................ .................. 204 6 3 Minority Carrier Mobility Fitting Parameters a t 300 K ................................ ....... 204 6 4 Overview of Simulation Variables ................................ ................................ ..... 205 6 5 Majority Carrier Mobility Fitting Parameters ................................ ...................... 205 6 6 Minority Carrier Mobility Fitting Parameters ................................ ...................... 205 6 7 Temperature Fitting Parameters. ................................ ................................ ...... 205 PAGE 10 10 LIST OF FIGURES Figure page 1 1 Single Event Effects terminology ................................ ................................ ........ 32 1 2 ................. 32 1 3 Simplified diagram of typical particle radiation spectra ................................ ....... 33 1 4 Particle composition of galactic cosmic rays. ................................ ...................... 33 1 5 Standard 6T SRAM in storage mode with a radiation event occurring ............... 34 1 6 Simulation results showing no ups e t (left) and upset (rig ht) for a 6T SRAM. ...... 34 1 7 Logic technology node and transistor gate length versus calendar year ............ 35 1 8 The problem of scaling. A lthough feature sizes are reduced, scaled devices are more susceptible to SEU ................................ ................................ .............. 35 1 9 T he variation of upset threshold with feature size for memory cells ................... 36 1 10 Observed MBU patterns from the testing of a 65 nm SRAM array. .................... 36 2 1 Illustration of the electron hole pair generation process. ................................ .... 59 2 2 Flowchart of the electron hole pair generation process. ................................ ..... 59 2 3 Radiation ionization energy ................................ ................................ ................ 60 2 4 Linear energy transfer vs. ion energy in silicon as calculated by SRIM .............. 61 2 5 Cylindrical Gaussian distribution for N ion ................................ ............................. 61 2 6 Electron hole density plot for a 590 nm single photon excitation process ......... 62 2 7 Electron hole pair distributions used in the simulations. ................................ ..... 62 2 8 Baseline simulation structure for the N+/P diode ................................ ................ 63 2 9 ......................... 63 2 10 Charge collection mechanisms of a particle strike in a reverse biased N+/P diode. ................................ ................................ ................................ .................. 64 2 11 effec t. ................................ ................................ ................................ .................. 64 PAGE 11 11 2 12 Typical shape of the single event charge collection current at a junction. .......... 65 2 13 Illustration of the depletion regi on width and the funneling depth. ...................... 65 2 14 Charge collection for the N+/P diode with the substrate doping varied. ............. 66 2 15 Potenti f ........................ 66 2 16 Charge collection results comparing diodes with an epitaxial configuration. ...... 67 2 17 Current transient results com paring diodes with an epitaxial configuration. ....... 67 2 18 Current transients for different LETs using the example N+/P Si diode. ............. 68 2 19 FLOODS predicted current transient for various constant mobility values. ......... 68 2 20 FLOODS predicted charge collection for various con stant mobility values. ........ 69 2 21 A reversed biased p n junction showing electron and hole currents ................... 69 2 22 Simulation results for a charge strike in a uniformly doped resistor with no bias applied. ................................ ................................ ................................ ....... 70 2 23 Simulation results for different carrier mobilities ................................ ................. 70 2 24 Illustration of the r ecombination process es ................................ ........................ 71 2 25 FLOODS predicted current transient with and without recombination. ............... 71 2 26 FLOODS predicted charge collection with and without recombination. .............. 72 2 27 Measured results of the Klaassen unified bandgap narrowing model ................. 72 2 28 FLOODS implemented model comparison of the Klaassen unified bandgap narrowing model ................................ ................................ ................................ 73 2 29 FLOODS predicted bandgap narrowing based on the Klaassen unified bandgap model. ................................ ................................ ................................ .. 73 2 30 Difference in current transients for the baseline diode simulation. ...................... 74 2 31 Difference in charge collection for t he baseline diode simulation. ...................... 74 3 1 Two dimension elements. ................................ ................................ ................... 86 3 2 Three dimension elements. ................................ ................................ ................ 86 3 3 Two dimension example for an area A i associated with a node ......................... 86 PAGE 12 12 3 4 NMOS I ON currents using the advanced physical model set for a variety of grid spacings. ................................ ................................ ................................ ..... 87 3 5 Percent change in solution time per Newton step for the FEQF when compared to the FVSG ................................ ................................ ...................... 87 3 6 An example of a perturbed mesh for the NMOS simulations. ............................. 88 3 7 The FVSG method loses accuracy for highly non Delaunay mesh conditions in the NMOS channel. ................................ ................................ ........................ 88 3 8 Single event upset comparison for both discretization methods. ........................ 89 3 9 Normalized average total transient simulation time. ................................ ........... 90 4 1 An increase in m grid points results in an increase in solution time. ................. 109 4 2 Flowchart of the proposed adaptive grid scheme. ................................ ............ 109 4 3 Example of grid refinement on a Gaussian function. ................................ ........ 110 4 4 Electron hole pair distributions used in simulations. ................................ ......... 110 4 5 N+/P diode 2 D simulation results comparing current transients ...................... 111 4 6 N+/P diode 2 D simulation results comparing co llected charge versus time ... 111 4 7 N+/P diode results. ................................ ................................ ........................... 112 4 8 Adaptive grid at peak refinement for the N+/P diode simulations. .................... 112 4 9 NMOSFET results. ................................ ................................ ........................... 113 4 10 Adaptive grid at peak refinement for the nMOS simulations. ............................ 113 4 11 Example of reflective symme try ................................ ................................ ........ 114 4 12 Illustration comparing: A) homogenous Neumann boundary. B) proposed diffusive boundary sink. ................................ ................................ .................... 114 4 13 Illustration showing the simulation structure with diffusive boundary sinks for two different widths. ................................ ................................ .......................... 115 4 14 Collected charge versus time for a reversed biased N+/P diode with reflective boundaries. ................................ ................................ ....................... 115 4 15 Collected charge versus time for a reversed biased N+/P diode with diffusive boundary sinks. ................................ ................................ ................................ 116 PAGE 13 13 4 16 Comparison of boundar y types with respect to device width and collected charge. ................................ ................................ ................................ ............. 116 4 17 Comparison of boundary types with respect to device width, collected charge and simulation time. ................................ ................................ .......................... 117 5 1 Ellipsoids of constant electron en .......................... 148 5 2 Simplified schematic of strain indu ced hole energy band splitting .................... 148 5 3 Linear spring element in equilibrium ................................ ................................ 149 5 4 Three dimensional stresses on an element. ................................ ..................... 149 5 5 Baseli ne tensor orientation notation and Miller indices for silicon. .................... 149 5 6 Piezoresistance factor P ( N,T ) as a function of impurity concentration and temperature for n type Si ................................ ................................ .................. 150 5 7 Laser induced current transient measurement system using a four point bending jig ................................ ................................ ................................ ........ 150 5 8 Schematic of N + /P diode structure through TEM/ EDS analys is and TEM image ................................ ................................ ................................ ................ 151 5 9 Schematic of laser induced current transients and 2 dimensional simul ation structure of an n+p diode ................................ ................................ .................. 152 5 10 Laser induced current transients as a function of uniaxial mechanical stress .. 152 5 11 Simulated laser induced current transients as a function o f uniaxial mechanical stress ................................ ................................ ............................. 153 5 12 Peak current as a function of mechanical stress. ................................ ............. 153 5 13 Co llected charges until 10 ns ................................ ................................ ............ 154 5 14 Orientation for the N+/P diode experiment and simulations. ............................. 154 5 15 Strained Si CMOS technology for 45 nm node ................................ ................. 155 5 16 Lattice expansion from germanium ................................ ................................ ... 155 5 17 TEM microg raphs of 45 nm node transistors ................................ .................... 155 5 18 2 D simulation structure. ................................ ................................ ................... 156 5 19 3 D MOSFET structure and Helium particle strike path ................................ .... 156 PAGE 14 14 5 20 Mea sured I V characteris tics for 45 nm strained Si CMOS .............................. 157 5 21 FLOODS predicted I D V GS characteristic for a strained silicon NMOS device .. 157 5 22 FLOODS predicted I D V DS characteristic comparing a strained and unstrained NMOS device ................................ ................................ ................................ ... 158 5 23 FLOODS predicted I D V GS characteristic for a strained silicon PMOS device .. 158 5 24 FLOODS predicted I D V DS characteristic comparing a strained and unstrained PMOS device ................................ ................................ ................................ .... 159 5 25 I D,lin enhancement versus un iaxial longitudinal tensile stress plotted for 10 and 0.1 ................................ ................................ .......................... 159 5 26 MOSFET orientation (and associated notation) with the channel in the <110> direction. ................................ ................................ ................................ ........... 160 5 27 2 D NMOS Stress XX component (c hannel direction) in [Pa] units .................. 160 5 28 2 D NMOS Stress ZZ component (depth direction) in [Pa] units ...................... 161 5 29 2 D NMOS current transient for strained and unstrained devices. .................... 161 5 30 2 D NMOS charge collection for strained and unstrained devi ces .................... 162 5 31 3 D NMOS Stress XX component (c hannel direction) in [Pa] units .................. 162 5 32 3 D NMOS Stress YY component (perpendi cular to channel) in [Pa] units ...... 163 5 33 3 D NMOS Stress ZZ component (depth direction) in [Pa] units ...................... 163 5 34 3 D NMOS current trans ient for strained and unstrained devices .................... 164 5 35 3 D NMOS charge collection for strained and unstrained devices .................... 164 5 36 2 D P MOS Stress XX component (channel direction) in [Pa] units. .................. 165 5 37 2 D PMOS Stress ZZ component (depth direction) in [Pa] units ....................... 165 5 38 2 D PMOS current transient for strained and unstrained devices .................... 166 5 39 2 D PMOS charge collection for strained and unstrained devices ................... 166 5 40 3 D PMOS Stress XX component (c hannel direction) in [Pa] units ................... 167 5 41 3 D PMOS Stress YY component (perpendicular to channel) in [Pa] units. ..... 167 5 42 3 D PMOS Stress ZZ component (depth direction) in [Pa] units ....................... 168 PAGE 15 15 5 43 3 D PMOS current transient for strained and unstrained devices. .................... 168 5 44 3 D PMOS charge collection for strained and unstrained devices .................... 169 5 45 Hysteresis effect of the deposited film as a function of temperature ................. 169 5 46 NMOS Stress XX component (channel direction) for STI induced stress. ........ 170 5 47 NMOS Stress YY component (perpendicular direction) for STI induced stress.. ................................ ................................ ................................ .............. 170 5 48 NMOS Stress ZZ component (depth direction) for STI induced stress. ............ 171 5 49 3 D NMOS current transient for STI strained and unstrained devices. ............. 171 5 50 3 D NMOS collected charge for STI strained and unstrained devices. ............. 172 5 51 Electron mobility change along the particle strike path in the <001> direction 172 5 52 Hole mobility change along the particle strike path in the <001> direct ion ....... 173 6 1 Qualitative comparison of commonly used bulk silicon mobility models for device simulation ................................ ................................ .............................. 206 6 2 Majority electron mobility as a function donor concentration for different mobility models at 300 K. ................................ ................................ ................. 206 6 3 Majority hole mobility as a function of acceptor concentration for different mobility models at 300 K. ................................ ................................ ................. 207 6 4 Majority electron mobility as a function of temperature and donor concentration. ................................ ................................ ................................ ... 207 6 5 Majority hole mobility as a function of te mperature and acceptor ..................... 208 6 6 Minority electron mobi lity in p type silicon at 300 K ................................ .......... 208 6 7 Minority hole mobilit y in n type silicon at 300 K. ................................ ............... 209 6 8 Sum of electron and hole mobility as a function of carrier concentration ve rsus experimental data at 300 K ................................ ................................ ... 209 6 9 Schematic of la ser induced current transients ................................ .................. 210 6 10 Electron hole pair distributions used in the simulations. ................................ ... 210 6 11 Simulated laser induced current transients in a reverse biased Si N+/P diode. ................................ ................................ ................................ ................ 211 PAGE 16 16 6 12 FLOODS predicted charge collection in a reverse biased Si N+/P diode. ........ 211 6 13 Simulated current transients in a reverse biased Si N+/P diode. ...................... 212 6 14 FLOODS predicted charge collection for a reverse biased Si N+/ P. ................ 212 6 15 Simulated current transients in a reverse biased Si N+/EPI/P+ diode. ............. 213 6 16 FLOODS predicted charge collection for a reverse biased Si N + /EPI/P + diode. ................................ ................................ ................................ ................ 213 6 17 Contributions to the majority electron mobility. ................................ ................. 214 6 18 Comparison of the propo electron mobility ................................ ................................ ................................ 214 6 19 for majority hole mobility ................................ ................................ ................................ ............. 215 6 20 Minority electron mobi lity in p type silicon at 300 K ................................ .......... 215 6 21 Minority hole mobi lity in n type silicon at 300 K ................................ ................ 216 6 22 Sum of electron and hole mobility as a function of carrier concentration ve rsus experimental data at 300 K ................................ ................................ ... 216 6 23 Comparison of electron mobility as a function of acceptor site and/or hole density. ................................ ................................ ................................ ............. 217 6 24 Majority electron mobility as a function of temperature and d onor concentration ................................ ................................ ................................ .... 217 6 25 Majority hole mobility as a function of temperature and acceptor concentration. ................................ ................................ ................................ ... 218 6 26 Schematic of the N/P/N simulation structure ................................ .................... 218 6 27 FLOODS 2 D simulation results for a saturation to cut off transi ent ................. 219 6 28 Schematic of laser induced current transients ................................ .................. 219 6 29 Simulated laser induced current transients in a reverse biased Si N+/P diode for the general purpose model ................................ ................................ ......... 220 6 30 FLOODS predicted charge collection in a reverse biased Si N+/P diode for the general purpose model. ................................ ................................ ............. 220 6 31 Enhanced Lombardi electron mobility model ................................ .................... 221 PAGE 17 17 6 32 Enhanced Lombar di hole mobility model overlaid on the measured hole mobility data ................................ ................................ ................................ ..... 222 6 33 Enhanced Lombardi e lectron mobility model overlaid on the measured electron mobility data ................................ ................................ ........................ 223 6 34 Electron and hole drift velocity in silicon as a function of electric field at three different temperatures. ................................ ................................ ..................... 224 PAGE 18 18 LIST OF ABBREVIATION S BJT Bipolar junction transistor CMOS Com plementary metal oxide semiconductor E H Electron hole EPI Epitaxial layer FEQF Finite element quasi Fermi FVSG Finite volume Scharfetter Gummel FLOODS Florida object oriented device simulator FLOOPS Florida object oriented process simulator ITRS Internati onal technology roadmap for semiconductors LET Linear energy transfer MBU Multiple bit upset MOSFET Metal oxide semiconductor field effect transistor NMOS N type metal oxide semiconductor field effect transistor PMOS P type metal oxide semiconductor field effect transistor PDE Partial differential equation SEE Single event effect SEU Single event upset SET Single event transient Si Silicon SPA Single photon absorption SRAM Static random access memory STI Shallow trench isolation TCAD Technology computer aid ed design PAGE 19 19 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy ENHANCEMENTS IN CMOS DEVICE SIMULATION FOR SINGLE EVENT EFF ECTS By Daniel J. Cummings December 2010 Chair: Mark E. Law Major: Electrical and Computer Engineering Single event effects in microelectronics can cause changes in memory state in spaceborne, airborne, and even terrestrial electronics due to the result ing charge collection from a radiation particle strike T he simulation of single event effects is an increasingly important area of numerical device simulation since the sensitivit y of microelectronics to single event upset is expected to increase as techn ology scaling continues. An especially important area of study for single event effects is in complementary metal oxide semiconductor (CMOS) transistor technology. As devices are downscaled a reduction in the amount of charge held on memory storage nodes increases CMOS vulnerabil ity to single event upset Single event upset experiment test costs are extremely high and require beam time at high energy ion accelerator facilities. Thus, d evice simulations are a useful way to predict and interpret device behav ior for such conditions, since comprehensive experimental testing for all particles, angles, and energies of interest is impractical. Many challenges exist in the area of single event device simulation. Firstly, modern technology computer aided design (TC AD) tools were not originally designed with single event simulations in mind. A particle strike generates a high density of PAGE 20 20 electron hole pairs along into the bulk of the device and often in non uniform patterns. Thus, gridding the simulation structure aro und the strike path requires significant TCAD expertise and the addition of grid points significantly increases solution time. Furthermore, the current flow around the strike path is isotropic in nature and is often not aligned with the device grid, making solution convergence problematic. Secondly, newer processing techniques such as strained silicon technology have continued to enable the scaling of CMOS devices by increasing carrier mobility. Process induced channel stressors such as embedded silicon ger manium and compressive and tensile capping layers introduce new complexities that need to be accounted for in single event simulations. Thirdly, the mobility models implemented in modern TCAD tools are inaccurate since they do not account for electron hol e scattering correctly. Because a high injection carrier condition occurs during a particle strike, the carrier scattering mechanism needs to be modeled accurately. This work addresses the challenges of single event simulation by presenting solutions to t he problems discussed above. First, a quasi Fermi finite element discretization approach is given to address the problems of single event simulation solution convergence and simulation time. Next, the problems associated with gridding around a particle str ike are discussed and an adaptive grid scheme is proposed. The proposed scheme offers a reduction in simulation time while retaining accuracy in results. Then, a piezoresistance mobility model is developed in order enable the single event simulation of str ained silicon CMOS devices. The results provide insight into the effects of strained silicon on charge collection. Finally, two new approaches to modeling electron and hole mobility are introduced to address the problem of electron hole PAGE 21 21 scattering in exist ing mobility models. Comparison tests show that the use of the new mobility models significantly improves the accuracy of the simulation results. The overall benefit of the above enhancements for the single event modeler is a savings in simulation time, an increased probability of solution convergence and an increase in accuracy. PAGE 22 22 CHAPTER 1 INTRODUCTION AND BAC KGROUND 1.1 Motivation Single event effects (SEE) in microelectronics occur when sensitive regions of a microelectronic circuit are struck by highly energetic pa rticles present in the natural space environment. For example, high energy heavy ions, alpha particles, protons, or secondary particles produced by neutron interactions can cause changes in memory state in spaceborne, airborne, and even terrestrial electro nics due to the resulting charge collection. T he simulation of single event effects is an increasingly important area of numerical device simulation since the sensitivit y of microelectronics to single event upset is expected to increase as technology scali ng continues [ 1 ]. An especially important area of study for single event effects is in complementary metal oxide semiconductor (CMOS) transistor technology. CMOS planar transistors have dominated the past two decades as the technology of choice for integra ted circuits (ICs) and a larger number of commercial ICs are being used in space and avionics applications. Advancements in process technology and a competitive electronics market have enabled 22 nm over the pa st 40 years [ 2 ] Co nsequently, as devices are downscaled a reduction in the amount of charge held on storage nodes increases device vulnerabil ity to single event upset [ 3 ] Single event upset experiment test costs are extremely high ( ~ $50 000 per part typ e ) and require beam time at high energy ion accelerator facilities [ 4 ] Thus, d evice simulations are a useful way to predict and interpret device behavior for such conditions, since comprehensive experimental testing for all particles, angles, and energies of interest is impractical. PAGE 23 23 1.2 Brief Overview of Single event Effects Radiation effects can have a large impact on the reliability of electronics in both the space and terrestrial radiation environments. Single events are name d as such because they depend o n the interaction of a single particle. This distinguishes them from other radiation effects (i.e. total ionizing dose) which depend on the dose or damage deposited by large number of particles. Single event effects can cause either rs or non single event effects that can cause malfunction in microelectronic devices. Figure 1 1 gives an overview of single event effects terms that are commonly used in indus try and Table 1 1 gives a description for each term. Most commonly, single event upset and latch up are the cause for malfunctions. The focus of this work is in the area of soft errors also known to as single event upsets (SEU) where a single particle stri ke causes a change in memory state. However, the simulation tool enhancements and physical model improvements presented in this work are also applicable and useful for all other soft error and hard error applications. The rate that soft errors occur is r eferred to as the soft error rate (SER) and the metric associated with SER and hard errors is referred to as failure in time (FIT). One FIT is equal to one failure per 10 9 device hours. For most electronic components the typical failure rate is about 20 20 0 FIT [ 5 ]. However, if mitigation and hardening techniques are ignored, the FIT can easily exceed 50,000 per chip. This can be very problematic for systems that require 100% uptime such as financial servers, commercial satellites and avionics equipment. As transistors are down s caled to meet consumer demand for faster, functional and efficient electronics, the device susceptibility to PAGE 24 24 radiation effects also increases dramatically. Therefore, it is important to understand the mechanisms for SEEs and device si mulation tools can be very useful in this regard. 1.2.1 Brief History of Single event Effects The first confirmed cosmic induced SEUs were reported in 1975 by Binder although the error levels were very low at that time [ 6 ]. As time continued, it became increasi ngly evident that that cosmic radiation was responsible for satellite subsystem soft errors and the first models for predicting soft error rates were formulated. In the late errors at ground level where the primary source o f radiation was found to be contaminated packaging materials [ 7 ]. The first reports of SEU from solar radiation sources such as protons and neutrons also began to be published. A high abundance of protons exist in the space environment making this discove ry of critical importance for the space electronics industry. to counter these problems, newer methods for hardening electronics were widely developed by industry [ 1 ]. During this period, much interested was generated in SEE due to critical errors caused by cosmic ions in the Voyager and Pioneer probes [ 8 ]. Additionally, with this knowledge, expensive retrofits were performed to mitigate SEEs for systems such as the Landsat D and Galileo systems [ 8 large number of commercial manufacturers began offering radiation hardened devices. However, with the increased use of commercial electronics in space and advancements in device technology came additional pro blems in maintaining system reliability. The downscaling of technology created new challenges for SEE since it was shown that scaling resulted in an increase in soft error susceptibility. An overview of spacecraft that have been impacted by SEE is given in Table 1 2. PAGE 25 25 In present day, device susceptibility to SEU continues to be a large issue as new developments such as strained Si CMOS, multi gate transistors, and SiGe based devices introduce new complexities in understanding SEU susceptibility. Additional ly, the rise in terrestrial soft errors in commercial electronics is becoming an increasing area of concern and has even become industry wide product reliability metric [ 4 ]. 1.2.2 Radiation Sources A general knowledge of the radiation environment is useful for u nderstanding the sources of radiation that cause single event effects. A common source for radiation fields. The external field results from the solar wind that is con tinually emitted by the sun and consists of plasma and ionized gas. The internal (or geomagnetic) field originates from within the Earth and is approximated by a dipole field. The trapped particles can be mapped in terms of the dipole coordinates that esti Charged particles are trapped by the magnetic field and then spiral and move along the magnetic field lines as in Figure 1 2 In addition to moving along the magnetic field lines, the trapped partic les drift longitudinally around the Earth where electrons drift eastward and protons move westward. The region is also known as the radiation belt environment [ 9 ]. Typical proton energies can reach several hundred MeV. Trapped protons are known to cause total ionizing dose (TID) effects, displacement damage (DD) effects, and single event effects. Electrons reach energies of a few MeV and contribute to TID effects, displacement damage effects, and charging/discharging effects. The electron charging/discha rging effects can be either spacecraft surface charging caused mainly by low energy electrons or deep dielectric charging caused by high energy electrons. PAGE 26 26 events (SPE) creat e large fluxes of energetic protons and other particles. SPE are unpredictable in time and occurrence, magnitude, and duration. These events are typically composed of solar protons and alpha particles, but can also include heavy ions, electrons, neutrons, and gamma particles. The composition and amount of particles for any given SPE varies greatly. Galactic cosmic rays (GCR) originate outside the solar system and have a highly variant particle energy spectrum. GCRs are believed to be remnants from supernov a explosions. Cosmic radiation includes heavy and highly energetic ions with energies in excess of 10 20 eV. Particles with such high energies have been detected on Earth and cause intense ionization along their tracks. In addition to the terrestrial and sp ace environment sources, radioactive contaminants in packaging materials can also be a source for SEE in microelectronics. A diagram is given in Figure 1 3 that shows the energies for particles such as trapped electrons, protons, alph as, and heavy ions [ 10 ]. Additionally, the particle composition of galactic cosmic rays is given in Figure 1 4 1.2.3 Example: Single event Upset in a 6T SRAM A SEU is a change of state caused by a radiation particle (e.g. heavy ions, alph a particles, protons, neutron s) that strik es a sensitive node in a micro electronic devi ce, such as those in a microprocessor or semiconduc tor memory If a strike occurs near a sensitive node of a circuit the resulting drift and diffusion carrier action wi ll create a large current and voltage transient spike. Both logic and memory (i.e. DRAM, SRAM) circuits in microelectronics are susceptible to single event effects. A simple way to illustrate a state change due to a particle strike is by using a six transi stor (6T) SRAM as an example. A standard 6T SRAM cell consists of two cross coupled inverters and two PAGE 27 27 word line enables, with a total of two PMOS and four NMOS transistors as in Figure 1 5 The cross coupling creates a regenerative feedback loop that maint ains the data state right and so on. Additionally, if a highly energetic particle strikes near the node storing ho le pairs is generated along the strike trajectory. Even though the e h pair cloud has a net charge of zero, the separation of these carriers due to high fields (i.e. funneling and depletion regions) results in a current transient at the node. Following the strike, the charge will collected causing a quick drop in the stored voltage on the left node as shown in Figure 1 2 If the PMOS on the left node cannot supply enough current to prevent the voltage on the left node with a state occurs. An example of the current and voltage change with respect to time for the SR AM upset is given in Figure 1 6 The SRAM example is just one of many possible SEEs that can occur due to a particle strike. Often, the state change due to a single bit upset will propagate through a logic circuit and cause a multiple bit upset (MBU). Addi tionally, a particle strike path with a low angle of incidence can traverse through multiple bit cells, causing an MBU. 1.3 CMOS Scaling and Susceptibility The study of single event effects in CMOS devices is incredibly important due to aw is the empirical observation that component density and performance of integrated circuits doubles every two years [ 11 ]. The downscaling of feature size (roughly analogous to CMOS gate length) is illustrated in Figure 1 7 Due to c ontinual advances in technology such as new processing techniques, device PAGE 28 28 However, scaling is a problem for SEE in microelectronics. As devices become get smaller and faster, they store less charge on critical circuit nodes. For example, a scaled MOSFET has a smaller volume and an ion strike that may not affect a large device will have a much large impact on a much smaller device as in Figure 1 8 The amount of charge (generated by particle strike) required to cause an upset is referred to as Q crit and will be discussed later. It has been shown that simple scaling rules predict an increase in s oft error susceptibility of about 40% per technology generation node [ 12 ]. An example of SEU susceptibility with respect to feature size is given in Figure 1 9 The problems of scaling extend to the circuit level. Multiple bit upsets (MBUs) are now more common due to the fact that it is more likely f or a single strike path to traverse many sensitive nodes. A recent study o bserved M B U patterns from the testing of a 65 nm SRAM array with a Kr ion (LET = 28.9 MeV cm 2 /mg), angled at 78.5 d egrees from normal, parallel to the n well [ 13 ]. The M B U patterns s how a constant string of upsets where the ion strike occurred shown by Figure 1 10 This shows that not only does scaling increase single event susceptibility at the device level, it also increases the chance of multiple upsets to occ ur at the circuit level. Although scaling limits are being approached for planar CMOS transistors, the $300 billion worldwide industry will be slow to change. It has been estimated that the time frame to implement a radically new device is roughly 30 year s. Additionally, s i licon CMOS technology is on course to offer a billion transistor chips for about $1 within the next decade which will be a very difficult price point to displace [ 2 ]. Thus, CMOS will continue to be the dominant form of nanotechnology fo r the foreseeable future. Since an PAGE 29 29 increasing number of spaceborne systems are using commercially available electronics suites which utilize CMOS technology, u nderstanding the impact of scaling (and associated processing techniques, materials, etc.) will b e key in SEE mitigation and hardening techniques for future spaceborne microelectronics. 1.4 Single Event Device Simulation Challenges Many challenges exist in the area of single event device simulation. Firstly, modern technology computer aided design (TCAD) tools were not originally designed with SEE simulations in mind. A particle strike generates a high density of electron hole pairs along into the bulk of the device and often in non uniform patterns. Thus, gridding the simulation structure around the stri ke path requires significant TCAD expertise and the addition of grid points significantly increases solution time. Additionally, the current flow around the strike path is isotropic in nature and is often not aligned with the device grid making solution co nvergence problematic. Secondly, newer processing techniques such as strained silicon technology have continued to enable the scaling of CMOS devices by increasing carrier mobility. Process induced channel stressors such as embedded silicon germanium (e S iGe) and compressive and tensile capping layers introduce new complexities that need to be accounted for in SEE simulations. Thirdly, the mobility models implemented in modern TCAD tools are inaccurate since they do not account for electron hole scatterin g correctly [ 14 ]. Because a high injection carrier condition occurs during a particle strike, the carrier scattering mechanism needs to be modeled accurately. This work addresses the challenges of SEE simulation by presenting solutions to each problem dis cussed above. First, a quasi Fermi finite element discretization approach is given to address the problems of SEU solution convergence and simulation PAGE 30 30 time. Next, the problems associated with gridding around a particle strike are discussed and an adaptive g rid scheme is proposed. The proposed scheme offers a reduction in simulation time while retaining accuracy in results. Then, a piezoresistance mobility model is developed in order enable the SEU simulation of strained silicon CMOS devices. The results prov ide insight into the effects of strained silicon on charge collection. Finally, two new approaches to modeling electron and hole mobility are introduced to address the problem of electron hole scattering in existing mobility models. Comparison tests show t hat the use of the new mobility models significantly improves the accuracy of the simulation results. The overall benefit of the above enhancements for the SEE modeler is a savings in simulation time, an increased probability of solution convergence and an increase in accuracy. 1.5 FLOODS Simulation Tool The simulation tool used for this work is the Florida Object Oriented Device Simulator (FLOODS) [ 1 5 ]. The presented algorithms, models and methods were implemented in FLOODS using the C++ and tcl/tk programmin g languages. FLOODS us es the d rift diffusion transport model and can use both finite volume and finite element discretization methods which will be described in detail in chapter 3 It also supports a variety of mesh element types for 2 D (triangular, rect angular) and 3 D (tetrahedra, bricks) simulations The corresponding process simulation tool called FLOOPS is used to simulate the process induced strained Si profiles and the n type/p type ion implantation distributions in Chapter 5. PAGE 31 31 Table 1 1. Single E vent Effects terminology description. Acronym Term Description SEU Single Event Upset Temporary change of memory or control bit SBU Single Bit Upset Single bit upset by one event MBU Multiple Bit Upset Several bits upset by the one event SEFI Single E vent Functional Interrupt Control path corrupted by an upset SELU Single Event Latch up Device latches in high current state SEGR/B Single Event Gate Rupture/Burnout Gate destroyed in MOSFET Table 1 2 Spacecraft for which single event effects hav e impacted [ 8 ]. Period Spacecraft 1970 1982 DE 1, Galileo, INSAT 1, intelsat IV, Landsat D, LES 8, LES 9, Pioneer Venus, SMM, Tires N, Voyager 1982 1990 AMTE/CCE, DSCS, ERBS, Galileo Lander, GEOS 6, GEOS 7, Geosat, GPS 9521, GPS9783, GPS9794, HUT, I US, MOS 1, OPEN, Shuttle, SPOT 1, TDRS 1, TDRS 4, UOSAT 2 1990 1997 ADEOS, COBE, ERS 1 (SEL), ETS V (SEL), EUVE, HST, HST STIS, Kitsat 1, NATO 3A, PoSAT 1, S80/T, SOHO, spot 2, SPOT 3, STS 61, Superbird, TDRS 5, TDRS 6, TDRS 7, Topex/Poseidon, UOSAT 2, UOSAT 3, UOSAT 5, WIND, Yahkoh BCS PAGE 32 32 Figure 1 1 Single Event Effects terminology Figure 1 2 Trapped particle behavior with respec agnetic field [ 9 ] PAGE 33 33 Figure 1 3 Simplified diagram of typical particle radiation spectra from the space environment. Figure 1 4 Particle composition of galactic cosmic rays Hydrogen (protons) and Helium (alpha particles) nuclei account for the vast majority of GCR flux where as h eavy ions comprise for only ~1% [ 16 ] PAGE 34 34 Figure 1 5 Standard 6T SRAM in stor age mode with a radiation event occurring on the lef t node near the NMOS drain [ 5] Figure 1 6 Simulation results showing no ups et (left) and upset (right) for a 6T SRAM [ 17 ] PAGE 35 35 Figure 1 7 Logic technology node and transistor gate length versus calendar year [ 2 ]. Figure 1 8 The problem of scaling. Although feature sizes are reduced, scaled devices are more susceptible to SEU since the mass and energy of ions stays constant. PAGE 36 36 Figure 1 9 The variation of upset threshold with feature size for memory cells [ 8 ]. As feature size decreases, the charge nee ded to create an upset decreases as well. Figure 1 10 Observed M B U patterns from the testing of a 65 nm SRAM array. Each cell is represented by a square. A single Kr ion causes direct charge collection a nd well collapse source injection for a large n umber of array cells (red) [ 13 ]. PAGE 37 37 CHAPTER 2 PHYSICAL MECHANISMS OF SINGLE EVENT EFFECTS 2.1 Introduction In order to describe the simulation tool and physical modeling enhancements in this work, it is essential to have an understanding of the physical mechanisms behind single event effects. This section describes these basic mechanisms, starting with the particle strike and carrier generation process. Subsequently, the charge collection transport mechanisms and related phys ics (recombination, mobility, bandgap narrowing) will be discussed. 2.2 Carrier Generation When a particle travels through a material such as silicon, it loses kinetic energy mainly through interactions with the lattice atoms and electrons of that material an d leaves a trail of ionization in its path. The incoming particles can be a heavy ions, protons or neutrons and usually have energies on the order of millions of electron volts. The energy from the incident particle is transferred into the material in the form of high energy electrons, photons and phonons. The process results in the ionization of electron hole pairs and is shown in Figure 2 1 and a flowchart is given in Figure 2 2 Two primary mechanisms contribut e to the stopping of a particle, electronic stopping (atomic electrons) and nuclear stopping (elastic scattering of lattice atoms). Electronic stopping is due to coulombic collisions between the incident ion and lattice electrons produce delta rays (a.k.a delta electrons). Delta rays are highly energetic electrons that scatter away from the original strike path. The subsequent lower energy collisions between the delta rays and crystal lattice atoms excite additional valence band electrons to higher energy bands since many empty states exist well PAGE 38 38 above the conduction band. The excited atomic electrons then thermalize energy by emitting photons and phonons of various energies. The high energy delta rays are e excited atomic electrons are primary between electrons, photons and phonons then cascades into lower and lower energies [ 18 ]. Nuclear stopping is the highly energetic (kin etic) displacement of lattice atoms, which in turn can lead to defects in the semiconductor lattice. The kinetic energy of the displaced atom is transferred to other lattice atoms and electrons which results in ionization. This cascading effect continues t o lower energies where energy transfers continue to ionize electrons and provide phonons to the lattice. Particles with a higher energy typically have a longer stop range in the target material than those of lower energy. After the nuclear and/or electron ic stopping of the incident particle, the semiconductor lattice begins to return to equilibrium. The electrons thermalize energy as they start to settle in the lowest available energy states in the crystal. What remains are electrons in the conduction band and holes in the valence band in equal pairs. The semiconductor is still neutral since both carriers have the opposite charge. It should be noted that the nuclear and electronic stopping mechanisms are a function of target material. It requires a differen t amount of energy to ionize an electron hole from material to material as shown in Figure 2 3 For example, SiC consumes much more energy per generated electron hole pair than silicon due to a wider bandgap. This means that for the s ame incident particle type and energy, the resulting generated electron hole pair PAGE 39 39 density will be much smaller for SiC than Si. This illustrates why the use of different materials in semiconductor processing has interesting implications for SEE mitigation and hardening. Device simulation tools start with the distribution of electron hole pairs from the strike and then simulate the movement using semiconductor physics transport models (e.g. drift diffusion, thermodynamic, hydrodynamic). It is important to n ote that the particle strike energy distribution needs to be in the form of electron hole pairs since device simulation tools do not simulation atomic level interactions. Particle physics tools such as GEANT4, NOVICE and MRED calculate the atomic level par ticle physics of a strike for a given species and target material using Monte Carlo methods [ 19 ]. They then output the electron hole pair distribution in a useable form for device simulation tools. The Monte Carlo approach involves the solving of the Boltz mann kinetic equation. Arguably, a TCAD hydrodynamic transport approach could be used to estimate the impact ionization, carrier temperatures and could more computationally efficient. However, to date, very little research has been done in this area becaus e the existing Monte Carlo tools are used as the standard. 2.3 Particle Strike Models Device simulation tools model particle strikes by approximating the electron hole pair distribution along the strike path with analytical models that are a function of partic le species, mass, energy and the target material type. As stated in the previous section, the device simulations start at a point where electron hole pairs are assumed to be thermalized. PAGE 40 40 2.3.1 Linear Energy Transfer The commonly used term in describing the ener gy loss of a particle (per unit length) in a material is called the linear energy transfer (LET). The LET is a function of The LET typically reported in units of MeV cm 2 /mg but can be converted into units of electron hole pairs per unit length using equation ( 1 1 ) as ( 1 1 ) i is the average electron hole ionization energy of the material in eV. In other i is the amount of energy required to create an electron hole pair in a material. For example, an LET of 1 MeV cm 2 /mg in silicon can be equ ated to 6.410 4 electron hole pairs per micrometer using equation ( 1 1 ) The parameters for i and the densities for various target materials are given in table 2 1. To calculate the total amount of charge Q that is ionized (in coulombs) during a strike, the LET in terms of electron hole pairs can be multiplied as ( 1 2 ) where q is the elementary charge of an electron (1.60210 19 C). Equation ( 1 2 ) can be implemented as a piecewise function for a particle that has an LET that varies with range. An example of the LET for various ions in silicon is given in Figure 2 4 where the data was taken from the Stopping and Range of Ions in Matter (SRIM) software package [ 20 ]. Typically the maximum LET for terrestrial events is below LET of ~13 MeV cm 2 /mg whereas LET for space events can be much higher. The sto pping range for a particular ion is a function of its energy and the target material. For instance, an Fe PAGE 41 41 ion with an energy of 1 MeV and 100 MeV will have an average stopping range of 0.86 m and 19.32 m in Si respectively. 2.3.2 Heavy Ion Modeling The modelin g of the electron hole distributions generated by heavy ions traversing through a material is a still a frequent area of discussion. For Monte Carlo simulators such as MRED, the initial ion track and resulting delta rays are modeled using a Gaussian approx imation of the dE/dx or LET values for a given ion species and material. Arguably, if one took the average e h pair distribution from thousands of Monte Carlo ion strike simulations, the results would start to take the form of a cylindrical Gaussian distri bution. Modern device simulation tools model heavy ion by using a temporal Gaussian that is a function of LET [ 21 ]. For most TCAD models, the carrier distribution of resulting from the ion strike is of the form ( 1 3 ) where N pk is the maximum peak carrier concentration in cm 3 r is the radial distance from the strike center r 0 and is the straggle. The term B(z) is a function of distance from the surface of the device. B(z) is typically a piecewise LET function used to model the variation of LET versus range effect for a particular ion. For example, the Bragg peak could be m odeled as a function of B(z) with information taken from SRIM. Frequently, a 50 nm 1/ e radius is used to determine the straggle for N ion since it represents an average lateral distribution for ion energies ranging from 1 to 100 MeV. An example electron hol e pair distribution for N ion is given in Figure 2 5 where the 1/ e radius is 50 nm and the LET is a constant 20 MeV /mg/cm 2 to a depth of 30 m. To date, to use PAGE 42 42 SRIM to estimate the LET, stopping range, and straggle for a given ion energy and target material. Alpha particle (Helium ion) modeling follows the same approach as equation ( 1 3 ) where SRIM can be used to estimate the characteristics for a specific energy. 2.3.3 Pulsed Laser Modeling The pulsed picosecond laser has become an important tool for use in single event effects experiments, especially in the area s ingle event upset and single event latchup [ 22 ]. Since heavy ions can be challenging to replicate in an experimental setting, pulsed lasers are frequently used to create conditions similar to those produced by an ion strike The pulsed laser technique exci tes the carriers in a semiconductor (via photons) using a tightly focused, above bandgap optical excitation [ 22 ]. Each absorbed photon generates a single electron hole pair. For a single photon absorption (SPA), the generated carrier density drops off expo nentially with distance from the target surface. Other techniques such as two photon absorption can inject carriers deeper into the substrate of a device. However, the pulsed lasers in the experiments performed later in this work use single photon absorpti on. Therefore SPA modeling will be the focus of this section. To model the electron hole pair distribution that results from a pulsed laser source, 22 ]. These equations define expressions for the laser beam irradiance as a function of depth in the semiconductor material. The radial dependence of the laser pulse irradiance is given by ( 1 4 ) PAGE 43 43 where N is the density of free carriers, P is the pulse power, and r is the distance from the center of the laser. The longitudinal dependence of the beam radius w(z) is defi ned as ( 1 5 ) where w o is the beam radius, z is the longitudinal (dep th) position relative to w o n is the linear index of refraction and is the wavelength of the light. With the heavy ion Gaussian model, the 1/ e radius is used as the radial distribution metric. In the case SPA, the common metric is the confocal parameter z 0 The z 0 param eter bounds the 1/ e contour and is defined as ( 1 6 ) where 2 z 0 defines the outer contour for which the beam is well col limated. Having defined the pulse irradiance I 0 and the longitudinal dependence of the beam radius w(z) the density of laser generated carriers in cm 3 as a function of depth can be defined by ( 1 7 ) ( 1 7 ) is given in Figure 2 6 where the electron hole density is shown for a 590 nm SPA process in Si with an energy of 4.2 pJ and a spot size diameter of 1.2 m. For the experimental work dis cussed in later sections, a cavity dumped dye laser with a wavelength of 590 nm, a pulse energy of 218 pJ, and a pulse width of 1 ps is used to inject electron hole p airs into a diode structure The laser direction is normally incident to the diode surface and ha s a spot size of 12 m in diameter. The electron PAGE 44 44 hole distribution generated by the laser in the experiments is shown in Figure 2 7 The carrier distribution for the experiments is more spread out due to the much larger spot si ze (12 m diameter ) than the example shown in Figure 2 6 (1.2 m diameter ). Table 2 2 gives the SPA model parameters that correlate to the experiment laser setup. 2.4 Charge Collection Mechanisms Charge collection is the primary mechanis m that causes single event effects. As discussed earlier, an SEE is more likely to occur if the energetic particle passes through a sensitive region of a microelectronic circuit. The previous section discussed the modeling of the electron hole pair distrib ution due to a particle strike. Once the generated electron hole pair distribution is known, the transport of these carriers can be solved with a device simulator. This section will discuss the physics behind charge collection. 2.4.1 Baseline Simulation Structu re To illustrate the mechanisms behind charge collection, a reverse biased N+/P diode structure is simulated for this section. A reversed bias N+/P diode is used because it represents the most sensitive regions of a modern microelectronic device (e.g. NMOS drain) and is more sensitive than a P/N diode [ 4] In a P/N diode, the charge collection is a function of hole mobility, which is much less than electron mobility. The 2 D simulation structure is 30 m by 40 m in width and depth. The N+/P diode accuratel y characterizes all the essential charge collection mechanisms (even in two dimensions) and is shown in Figure 2 8 To mimic an ion strike, the electron hole distribution is modeled using equation ( 1 3 ) and has a constant LET of 1 MeV cm 2 /mg. The peak carrier concentration of the strike is 8.2110 18 cm 3 has a 1/ e radius of 50 nm PAGE 45 45 and terminates at a depth of 20 m. The doping profile is g iven in Figure 2 9 for the 2.4.2 The Basics of Charge Transport The three physical mechanisms that determine charge transport after a strike are the drift, diffusion, and recombination of carriers as shown in Figur e 2 10 The figure shows the worst case path for a strike since it traverses the depletion region where a high field exists. At the beginning of charge collection process, a cylindrical track of electron hole pairs at a very high conce ntration is formed along the strike path. The high field inside the depletion region of the reversed biased device is very effective at collecting the charge through the drift process. Prior to the particle strike the majority of the voltage drop exists a cross the depletion region. The high injection of electron hole pairs temporarily eliminate s the depletion region and most of the voltage drop occurs over the area in the vicinity of the ion track. In other words, the high injection carrier distribution a long the strike path will extend the junction field deep into the device since the highly conductive charge neutral plasma is high enough in density to disturb the local field The high field in the previous depletion region redistributes around the vicini ty and bottom of the strike track in the form of a funnel. A good analogy would be to think of the depletion region extending down and around the strike path (temporarily). Thus, carriers that are created further from the original depletion region in the i nitial strike track will be collected This disturbance to the junction electrostatic potent ial is known as the funneling The c arriers in the track remain in a vertical field and separate. For the case of the reverse biased N+/P diode, e lectrons drift up to the positive potential and holes PAGE 46 46 drift down to the substrate [ 18 ] The funneling effect can be seen by the migration of electrostatic potential contours as shown in Figure 2 11 for the baseline simulation. Following the drift acti on and the collapse of the funnel, the remaining carriers continue to diffuse where they are then collected in the depletion region or substrate via contacts. In addition to the drift/diffusion transport, the number of carriers is reduced over time by reco mbination. The typical current transient shape for the funnel creation, drift and diffusion transport mechanisms is shown by Figure 2 12 2.4.3 Analytic Approximations Before TCAD tools were widely available, analytic equations were used t o predict SEE behavior for devices. A common predict ion for the depth d of fu nn el collection in an N +/ P junction below the N+/P junction edge is given as ( 1 8 ) where n,p are the electron and hole mobilities and W is the depletion region width after the funneling effect ends as in Figu re 2 13 If one assumes that the electron mobility is twice the value of the hole mobility, equation ( 1 8 ) reduces to the funnel depth being equal to three depletion layer wi dths. To approximate the shape of the single event current pulse (shown in Figure 2 12 ), Messenger developed a model for the pulse in the form of a double exponential given by ( 1 9 ) where N is the electron hole pairs per unit length, E 0 the maximum field, the high l ) is t he time constant of charge collection from the funnel and PAGE 47 47 (sec l ) is the time constant for the initial formation of the funnel region [ 23 ]. This formulation is often implemented in circuit simulations since it is in a friendly form to be used as a curren t source. An analytical model of the funnel effect on total collected charge was developed by McLean and Oldham [ 24 ]. The model assumes that the temporal and spatial history of the funnel field can be estimated using an effective field that is related to the relaxed depletion region field after the event. This leads to the following equation for collected charge ( 1 10 ) with the collection time as ( 1 11 ) where D is the diffusion constant, v p is the escape velocity for holes, N 0 ,avg is the average carrier density along the track, and N 0 is the density near the surface. Although this model overpredicts collected charge, it is useful for firs t order estimations. For example, if the mobility n is increased, more charge is collected. 2.4.4 Doping Profile Effects The doping profile of a structure has a direct impact on charge collection for several reasons. First, it determines size and duration of t particle strike. Also, the doping profile (and external bias) determines the size of the depletion region. A larger depletion region cross section means it will be more likely for a carrier to diffuse into the region and then be collected at a contact. PAGE 48 48 Previous work has shown that for an N+/P diode structure, a lighter substrate doping results in a longer funneling depth and thus, an increase in collected charge [ 15 ]. For example, Figure 2 14 shows the impact of substrate doping for a reverse biased N+/P diode. The only difference in this example, with respect to the baseline simulation structure described earlier, is that the entire p type region is uniformly doped (no p well). As shown by the figure, a higher substrate doping results is less collected charge. This is because the high density charge region will perturb the local field less, due to the higher background doping. In other words, for a strike in a highly doped region, the non equilibrium ch arge density is less versus the steady state density ( n,p >> n 0 p 0 ) than for the lightly doped case. For the lightly doped case, the funnel collects most of the generated charge via drift as in Figure 2 15 This is because most o f the generated charge (strike path 20 m deep) falls inside the funnel region which leaves little to be collected by diffusion. For every case, as charge separates out, the funnel begins to collapse. However, the funnel lasts longer for the light doped re gion since it takes more time for the generated carrier density to fall to the pre strike density. This characteristic has been discussed in detail by Hsieh [ 25 ]. A common method to implement the advantage of higher doping levels for SEE is to use a highl y doped epitaxial layer in the substrate. For an N+/P diode, an EPI structure with a heavily doped p type substrate will limit the funneling length to the depletion region between the EPI and N+ region [ 15 ]. For this comparison, a N+/P sub/P+ EPI diode str ucture with doping levels of 10 20 cm 3 (N+), 10 16 cm 3 (P substrate), and 10 18 cm 3 (P+) was created. The junction depths for the EPI diode were 0.2 m (N+/P sub) and 2.5 m (P sub/P+) respectively. Results comparing the EPI diode and PAGE 49 49 an N+/P diode (p sub = 10 16 cm 3 ) are shown in Figure 2 17 and Figure 2 18 Since the highly doped EPI hinders funnel formation, the charge collection due to drift is limited to the lightly doped p type region. This results i n less charge being collected for the EPI structure. This is common technique for mitigating charge collection (and single event latch up) in CMOS devices at the expense of a more complicated front end process. Although less total charge is collected, the magnitude of the peak current is for the EPI diode higher for the initial portion of the single event transient These examples show the impact of doping profiles on single event behavior. However, many other doping methods exist and can be implemented for single event hardening. For example, retrograde wells and buried layers can be used to create internal electric field s that change how charge is collected. The energy of a particle is proportional to the linear energy transfer and subsequently, the amoun t of electron hole pairs generated along the particle strike path. For example, if the LET value for the cylindrical Gaussian equation ( 1 3 ) is increased, an increase in collected charge is observed. The impact of various laser energies on the current transient is shown in Figure 2 18 2.4.5 Mobility The results of semiconductor device simulations are highly dependent on the electron and hole mobility models. For instance, the overall effect of mobility on current density can be shown in terms of quasi Fermi levels as ( 1 12 ) ( 1 13 ) PAGE 50 50 where n and p are the electron and hole densities, n,p the quasi Fermi levels, J n,p the n,p the mobilities. Therefore, it is important to choose an accurate mobility model so that the simulation results will be relevant. Advance d mobility models will be discussed in great detail in chapters 5 and 6. However, the baseline example diode uses the assumption that electron mobility is a constant 200 cm 2 /Vs and hole mobility is a constant 100 cm 2 /Vs. To show the impact of mobility on current transients and charge collection, the baseline values for mobility were multiplied by factors of one half, two, and three. Figure 2 19 and Figure 2 20 show the impact where it can be seen that the mobili ty constants are proportional to the amount of collected charge. Interestingly, Figure 2 20 shows a difference in collected charge even though the particle strike LET value was constant. Thus, w hen examining the above results, th e question arises as to why there is a difference in collected charge if the e h pair distribution is the same for each simulation. The difference is due mainly to the funneling mechanism. During the funneling process, the mobility value dictates how fast (and thus how much) charge is swept to the contact. The mobility also impacts how many carriers are swept into the depletion region via diffusion, though this is a secondary effect. The next section will discuss in great detail how charge is conserved in a device during a single event. 2.4.6 Charge Conservation To explain charge conservation in a device during a single event, take for example Figure 2 21 where a reverse biased N+/P junction is shown. If an electron hole pair is created in t he p type material within a hole diffusion length of the depletion region, the hole may diffuse left, get caught in the drift field and pass through the p type region without recombining. At the far left contact, the hole would then recombine with an PAGE 51 51 elect ron pulled off the wire. Additionally, if an electron and hole arrive at an ohmic contact at the same time, they are annihilated by recombination. If two electrons and one hole arrive at the contact, only one electron would be collected. Take for example a charge strike in an unbiased, uniformly doped resistor. Assume ohmic contacts are placed on the left and right bounds of the resistor and that the carrier mobility is a simple constant. A particle strike the middle of the resistor would then generate a la rge number of electron hole pairs. Since there is no applied field and the mobility is constant, the electrons and holes would diffuse at the same rate. The electrons and concen tration. Therefore, the net current collected is zero since the electron and hole flux is always the same at the contacts as shown by the simulation results in Figure 2 22 This idea can be extended to the simulation results of the p revious section. If we look at the mobility results for the one half mobility factor and the three times mobility factor in Figure 2 20 we see a large difference in collected charge. However, if we sum the total electron and hole cur rent in both the top and bottom contacts (and neglect recombination at the contacts and in the device) the same result in collected charge is observed as shown in Figure 2 23 In fact, the same amount of charge is collected at the con tact as is deposited in the device initially. However, due to the factors such as mobility and doping levels, only a portion of the deposited charge is collected at the top contact. Therefore, for a given amount of charge generated, the mobility, funnel de pth (doping), funneling time and other parameters can influence the ratio of charge PAGE 52 52 collected at the N+ junction. This is why a difference in collected charge is observed when using different mobility values, even if the LET is constant. 2.4.7 Recombination A l arge number of electron hole pairs are generated during a particle but not every e h pair will reach a contact due to the possibility of recombination. Recombination is a built in characteristic of semiconductor devices that acts to reduce the charge colle ction. When a semiconductor is perturbed from a state of equilibrium, it has an excess or deficit of carriers relative to their equilibrium values. Recombination generation (R G) acts as the order restoring mechanism that seeks to stabilize or eliminate th e perturbation [ 26 ]. Since non equilibrium conditions exist during normal device operation, recombination generation will always have an influence on device characteristics. For the case of single event effects, carriers from a particle strike create an ex cess of carriers relative to the equilibrium state. Therefore, recombination is an important mechanism to model for SEE simulations. As the name implies, when an electron and hole are pulled together by coulombic forces, the conduction band electron can e nter the empty valence band state and recombine. The recombination event conserves energy such that if an electron recombines, energy must be released in the form of photons or phonons. The recombination rate varies between high level and low level injecti on levels. The primary mechanisms in silicon are Shockley Read Hall (SRH) R G center recombination and Auger band to band recombination. Although other recombination mechanisms may exist, their effects are considered insignificant for silicon although furt her studies would be beneficial [ 18 ]. PAGE 53 53 2.4.7.1 Auger Recombination In the Auger process, band to band recombination occurs when two like carriers collide. The energy released by the recombination mechanism is transferred to the remaining carrier as in Figure 2 24 energy and the other electron recombines. The equation that defines Auger recombination is given by ( 1 14 ) where C n,p are temperature independent coefficients. The temperature dependent coefficients can be written as ( 1 15 ) where the subscripts i ( n,p ) stand for electrons or holes. The standard coefficient values are listed in T able 2 process thorough a potential barrier (i.e. Zener process). 2.4.7.2 SRH Recombination SRM recombination is the transition of electrons and/or holes to states (R G centers) near the middle of the bandgap. Common impurities with near midgap energy levels are Au, Cu, Mn, Cr, and Fe. The recombination at an R G center is a two step process. For example, a hole could come into the vicinity of a trapped electron, become attracted to the electron, lo se energy and then annihilate the electron within the center. Alternatively, an electron can lose energy a second time from a midgap state and annihilate a valence band hole as shown in Figure 2 24 It should be noted that the R G PAGE 54 54 ce nter process is not limited only to near midgap energy states. SRH recombination is formulated as ( 1 16 ) with ( 1 17 ) and ( 1 18 ) where E trap represents of the difference in bandgap (eV) between the defect and intrinsic levels. The doping dependence of the model is ( 1 19 ) where the subscript i stands for e (electrons) or h (holes). The t ypical values for N ref n,max and p,max are given in Table 2 4. 2.4.7.3 The Impact of Recombination on Charge Collection Recombination plays an important role in charge collection. Using the bas eline device described earlier, the simulation results with and without recombination (SRH and Auger) are shown in Figure 2 25 and Figure 2 26 The results show that if recombination is neglected, the error in co llected charge is ~18%. This shows the importance of accounting for recombination effects in SEEs. PAGE 55 55 2.4.8 Bandgap Narrowing The energy bandgap in silicon E g narrows as a function of impurity concentration. This is due to the fact that concentration at high impur ity concentrations the density of energy states no longer has a parabolic energy distribution and becomes dependent on the impurity concentration [ 27 ]. This can have implications for single event behavior since particle strike paths often traverse highly d oped regions. Bandgap narrowing models for both n type and p type materials were developed separately by Slotboom and del Alamo [ 27 ],[ 28 ]. Subsequently, Klaassen formulated a unified bandgap narrowing model that works for both n and p type regions using only one set of model parameters [ 29 ]. The Klaassen bandgap model is implemented in FLOODS for this work. In this model, the effective intrinsic carrier concentration is given by ( 1 20 ) where T is temperature and g is the apparent bandgap narrowing. The bandgap narrowing term is independent of temperature and is defined by ( 1 21 ) where N is the impurity concentration and the remaining terms are fitting parameters. The parameters for the Slotboom, del Alamo, and Klaassen bandgap models are given in Table 2 5. A comparison of the expe rimental data and the models is given in Figure 2 27 and Figure 2 28 This shows that the Klaassen model results fall between the Slotboom and del Alamo approaches [ 29 ]. Figure 2 29 shows the predicted bandgap narrowing for the Klaassen model for the baseline example diode whereas impurity PAGE 56 56 concentration increases, the bandgap decreases. Figure 2 30 and Figure 2 31 show simulation results using the baseline example diode, with and without bandgap narrowing active. As shown by the figures, bandgap narrowing can play an important role in estimating collected charge and current transients for single event effects. 2.5 Summary This chapte r gave detailed descriptions of the physical mechanisms behind single events. Significant error can be introduced into the simulation results if any of the mechanisms are incorrectly modeled. Starting with the electron hole pair generation t he physics of carrier ionization and thermalization were described and equations that model particle strike carrier generation were discussed. The physics behind charge collection mechanisms such as drift, diffusion and funneling were explained and analytic equations fo r estimating the total charge collection and current transients were given. Next, the effects of doping, particle energy, mobility, recombination and bandgap narrowing on single event effects were discussed. All of the above mechanisms play an important ro le in single event effects and should be accurately modeled in modern device simulation tools. PAGE 57 57 Table 2 1. Parameters for calculating LET for various semiconductor materials [ 18 ]. Target Material eV / e h pair Density (gm/cm 3 ) fC/ m for an LET=1 MeV/mg/cm 2 Si 3.6 2.32 10.4 GaAs 4.8 5.32 17.8 InP 4.5 4.81 17.1 In 0.47 Ga 0.53 As 2.9 5.49 30.3 SiC 8.7 3.21 5.9 Table 2 2 Experimental parameters for the single photon absorption pulsed laser for silicon. Parameter Value Unit(s) Notes c 2.998e8 m/s Speed of light 1.0546e 34 Js 0.42e 12 s Pulse time 5.9e 7 m Wavelength 0.5824112 1/ m Absorption E 13.5e 12 J Laser energy w 0 6.0e 6 m Spot size diameter P W Power Table 2 3 Standard coefficients for Auger recombination model [ 21]. Parameter A A [cm 6 s 1 ] B A [cm 6 s 1 ] C A [cm 6 s 1 ] H [1] N 0 [cm 3 ] Electrons 6.7e 32 2.45e 31 2.2e 32 3.46667 1e18 Holes 7.2e 32 4.5e 33 2.63e 32 8.25688 1e18 PAGE 58 58 Table 2 4 Standard coefficients for SRH recombination model [21]. Parameter N ref [cm 3 ] max [s] Electrons 110 16 110 5 Holes 110 16 310 6 Table 2 5 Parameters for silicon bandgap narrowing models. Parameter Slotboom [27] (p type) del Alamo [28] (n type) Klaassen [29] (n and p type) C 1 (cm 6 K 3 ) 9.6110 32 1.3810 33 9.6110 3 2 C 2 0 .5 0 0.5 E g (V) 1.206 1.206 1.206 V 1 (mV) 9.0 9.35 6.92 N 2 (cm 3 ) 1.010 17 7.010 17 1.310 17 PAGE 59 59 Figure 2 1 Illustration of the electron hole pair generation process. Figure 2 2 Flowchart of the electron hole pair generation process [30] PAGE 60 60 Figure 2 3 Radiation ene rgy (from photons, hot elec trons, particles) consumed per generated electron hole pair, as a functi on of the bandgap width E G [30 ] PAGE 61 61 Figure 2 4 Linear energy transfer vs. ion energy in silicon as calculated by SRIM [ 20 ]. Note the blue region is the LET range for terrestrial events. Figure 2 5 Cylindrical Gaussian distribution for N ion LET = 20 MeV cm 2 /mg, 1/e radius = 50 nm. PAGE 62 62 Figure 2 6 Electron hole density plot for a 590 nm single photon excitation process in silicon for a 4.2 pJ, 1 ps pulse focused to a spot diameter of 1.2 m. The carrier density is plotted in electron hole pairs/cm3 [ 22 ]. Figure 2 7 Electron hole pair distributions used in the simulations. Single photon absorption, laser energy = 13.5 pJ, radius = 6 m PAGE 63 63 Figure 2 8 Baseline simulation structure for the N +/P diode used in this section. The ion strike path is directly in the center of the structure where the grid is dense. Figure 2 9 PAGE 64 64 Figur e 2 10 Charge collection mechanisms of a particle strike in a reverse biased N+/P diode [ 5] A B Figure 2 11 FLOODS predicted potential contour def effect. A) Profile near beginning of strike. B) Profile as time progresses. Strike path PAGE 65 65 Figure 2 12 Typical shape of the single event charge collection current at a junction. Figu re 2 13 Illustration of the depletion region width W and the funneling depth d. PAGE 66 66 Figure 2 14 Charge collection for the N+/P diode with the substrat e dopin g varied A B Figure 2 15 substrate doped at 10 16 cm 3 B) P substrate doped at 210 1 7 cm 3 PAGE 67 67 Figure 2 16 Charge collection results comparing diodes with an epitaxial N+/P sub/P+ and a N+/P sub (10 16 cm 3 ) configuration Figure 2 17 Current transient results comparing dio des with an epitaxial N+/P sub/P+ and a N+/P sub (P sub=10 16 cm 3 ) configuration PAGE 68 68 Figure 2 18 Current transients for different LETs using the example N+/P Si diode Simulations used a c yli ndrical Gaussia n distribution to generate e h pair profile Figure 2 19 FLOODS predicted current transient for various constant mobility values. PAGE 69 69 Figure 2 20 FLOODS p redicted charge collection for various constant mobility values. Figure 2 21 A reversed biased p n junction showing electron and hole currents in semiconductor and electr on currents in the circuit [ 26 ]. PAGE 70 70 Figure 2 22 Simulation results for a charge strike in a uniformly doped resistor with no bias applied. Figure 2 23 Simulation results for different c arrier mobilities showing the sum of collected electron and hole charge at the diode contacts (recombination neglected). Collected charge equals the deposited charge. PAGE 71 71 A B Figure 2 24 Illustration of the A) SRH recombination process and B) Auger band to band recombination process [ 26 ]. Figure 2 25 FLOODS predicted current transient with and without recombination PAGE 72 72 Figure 2 26 FLOODS predicted charge collection with and without recombination Figure 2 27 Measured results of the Klaassen unified bandgap narrowing model versus the Slotboom and de l Alamo models. [ 29 ] PAGE 73 73 Figure 2 28 FLOODS implemented model comparison of the Klaassen unified bandgap narrowing model versus the Slotboom and del Alamo models. Figure 2 29 FLOODS predicted bandgap narrowing based on the Klaassen unified bandgap model. PAGE 74 74 Figure 2 30 Difference in current transients for the baseline diode simulation. Results shown with (b aseline) and without bandgap narrowing effects. Figure 2 31 Difference in charge collection for the baseline diode simulation. Results shown with (baseline) and without bandgap narrowing effects. PAGE 75 75 CHAPTER 3 DISC RETIZATION METHODS F OR SEE SIMULATIONS 3.1 Introduction Three coupled nonlinear partial differential equations (PDEs) form the foundation of modern semiconductor device simulation [ 31 ]. These equations, consisting of the Poisson, electron continuity and hole continuity equations, can be solved using a variety of approaches [ 32 ]. The earliest work in device simulation started with Gummel, who simulated one dimension bipolar transistors by sequentially solving the three PDEs using the drift diffusion transport e quations an d an iterative solution method Building on this work, Scharfetter and Gummel demonstrated the stable upwind discretization of the transport equations [ 33 ]. This remains the most commonly used method (a.k.a. Scharfetter Gummel method) in modern device simulation tools since it is relatively computationally efficient, well tested and accurate. Most device simulators solve the three PDE equations by using electron density, hole density, and electrostatic potential as the solution variables ( n p ) and typically use a finite volume Scharfetter Gummel (FVSG) discretization scheme. However, discretization methods are not limited to just these variables. An alternate approach to the FVSG scheme is to solve the PDEs in terms of electron and hole quasi Fermi levels and electrostatic potential ( n p ) using a finite element approach [ 34 ],[ 35 ]. As will be shown in the simulation results, this finite element quasi Fermi (FEQF) approach holds several advantages over the FVSG approach for single event si mulations [ 31 ]. In a FVSG scheme, the current flow is computed on the grid edges. For most semiconductor devices, the simulation grid will already be aligned in the direction of current flow (e.g. MOSFET channel). However, for single events, a particle st rike PAGE 76 76 generates carriers throughout the device and the resulting funneling, drift and diffusion current is rarely aligned with the grid. For these non ideal conditions, the FEQF method could be more accurate and stable than the FVSG approach, since current flow in the FEQF method is not defined on the grid edges. Thus, it is important to compare the FVSG method to the less prevalent finite element quasi Fermi (FEQF) approach for single event simulations. The following section will start by describing the FVS G and FEQF discretization methods in detail. Next, the grid (or mesh) element types and the physical models implemented in the simulations will be defined. Then, simulation results from an NMOS device are used to show that both discretization methods give comparable results for an ideal grid and different results for a perturbed grid. Finally, a set of particle strike simulations are used to show the benefits of the FEQF method as compared to the FVSG approach in terms of convergence and simulation time. 3.2 Di scretization Overview The set of coupled, time dependent partial differential equations that govern semiconductor device behavior can be written as ( 1 22 ) ( 1 23 ) ( 1 24 ) where is the dielectric permittivity, q the elementary charge, is electrostatic potential, n and p are the electron and hole densities, N D + and N A are the ionized donor and acceptor densities, J n and J p are the electron and hole current densities, and U p a nd U n PAGE 77 77 are the net electron and hole recombination rates. To obtain a closed system of equations, the current densities are written as quasi linear functions of driving potential in gradient form ( 1 25 ) ( 1 26 ) where n p are the electron and hole quasi Fermi levels and n p are the mobilities. The quasi Fermi levels are functions of the electrostatic potential and the electron and hole carrier de nsities. For example, in the case of a nondegenerate semiconductor, the quasi Dirac) as ( 1 27 ) ( 1 28 ) where kT/q is the thermal voltage and n i is the intrinsic carrier concentration. Using these relations, the current density in equations ( 1 25 ) and ( 1 26 ) can be rewritten in the familiar relationship as the sum of drift and diffusion components ( 1 29 ) ( 1 30 ) where E is the electric field and D n,p is the diffusion coefficient Because the system has three PDEs and only three solution variables ( n p ), numerical approaches have to be taken to find the solutions in both time and space. These numerical approaches involve the discretization of the problem domain on a set PAGE 78 78 of predetermined, discrete points known as the grid (or mesh) using a set of algebr aic relations derived from equations ( 1 22 ) through ( 1 24 ) Since discretization methods are dependent on the geometry of the grid a variety of grids composed of different element types are used for two dimensional ( 2 D ) and three dimensional ( 3 D ) simulations. For 2 D simulations, triangular rectangular (quad) a nd pentagonal (five point) elements are typically used as shown in Figure 3 1 The five point element is used for terminating lines on a grid. For the generation of 3 D grids tetrahedron hexahedron (brick) and prism element types ar e used as in Figure 3 2 It should be noted that 3 D grids present significant challenges for discretization especially in the coupling of equations. For instance, to convert a hexahedral elements into tetrahedral requires dividing op posite faces on the hexahedra with different diagonals and then pulling tetrahedral out of the four corners (thus changing the coupling and bandwidth). Hexahedra and tetrahedral will be compared in this chapter since they are commonly used for 3 D simulati ons. Prisms are ignored for this work since they are only suitable for problems which have weak three dimensional effects [ 36 ]. 3.2.1 Finite Volume Discretization The finite volume method is used for representing and evaluating partial differential equations in the form of algebraic equations whose values are calculated at on a grid. When using electron and hole densities as solution variables in a finite volume scheme, each partial differential equation is integrated over a control volume A i surrounding each node as in Figure 3 3 The control volume C i is defined by the perpendicular bisec tors of the grid element sides. The Poisson, electron continuity and PAGE 79 79 hole continuity equations are written using divergence operators and are in the general form of ( 1 31 ) where F(x,y) is some vector function and u(x.y) is some scalar function. The divergence can be discretized on the grid as ( 1 32 ) where A i is the volume associate with node i as defined by the bounding line C i a nd n is an outward unit vector normal to C i explicitly guarantees conservation of carriers and charge. For example, the Poisson equation can be written using the form of equation ( 1 32 ) as ( 1 33 ) The evaluation of the electric field can be done as a straight line approximation across the edge, simplifying the process. The current J n,p is then be evaluated using the Scharfetter Gummel formula [ 33 ]. 3.2.2 Finite Element Discretization For t he finite element quasi Fermi scheme, the continuity equations can be rewritten in terms of equations ( 1 25 ) and ( 1 26 ) and the gradient of n,p can be computed over each mesh element. This means that the solution variables are now the quasi Fermi levels and electrostatic potential ( n p ). In finite element methods, the variational form of the problem is derived. For example, the variatio nal form of the Poisson equation in 2 D dictates that must satisfy the condition PAGE 80 80 ( 1 34 ) with the associated boundary conditions for all v(x,y) are ( 1 35 ) ( 1 36 ) ces of the domain which satisfy equation ( 1 34 ) [ 36 ]. In FLOODS, the discretization for the finite element method starts er subspaces (e.g. triangles or rectangles). Then, the subspaces are discretized into a set of points on which piecewise linear polynomial interpolation is used. In summary, this method is a process for producing an optimal piecewise interpolant to the tru e solution. As with the FVSG method, the FEQF is a function of grid points and grid density. 3.3 Simulation Methodology A variety of physical models are implemented in the simulation tool. For mobility, the simulation models used are the Philips mobility, Lom bardi mobility, and velocity saturation models. The mobility models are described in great detail in chapter 6. The Philips unified mobility model unifies the description of majority and minority carrier bulk mobilities and takes into account carrier carri er scattering, screening of ionized impurities, and clustering of impurities [ 37 ]. The Lombardi model is a function of surface acoustic phonon scattering and surface roughness scattering [ 38 ]. For recombination generation the Shockley Read Hall (SRH) and A uger band to band models were used. The physical models were divided into two sets for the simulations so that the discretization methods could be tested under different circumstances Table 3 1. In PAGE 81 81 addition, this comparison will show the impact of electric field dependent models such as Lombardi mobility and velocity saturation. 3.4 Simulation Results Both FVSG and FEQF methods were compared for a variety of mesh element types and device structures. Additionally, the x y z axis grid spacings were varied sin ce smaller spacings give a more accurate result but require more computation time. The assembly time, linear solve time, and number of Newton iteration steps were measured for each simulation. 3.4.1 Short Channel MOSFET results For a baseline comparison, a moder n N MOS device was simulated to compare the FVSG and FEQF methods. A simple short channel N MOS device with a 40nm gate and Gaussian doping profiles in the source/drain was created as a template. For 2 D N MOS simulations, the FVSG and FEQF methods performed very similarly in terms of output current and number of Newton steps required for convergence. For rectangular and quad diagonal mesh elements, both methods gave similar nMOS ION currents at very tight grid spacings ( Figure 3 4 ). The currents began to vary as the grid spacings increased above 1 nm, though both methods followed a similar trend. The assembly time for the FEQF approach was longer than the FVSG, and on average resulted in a ~22% increase in total solution time per Newton step ( Figure 3 5 ). The time increase is due to the fact that in the FVSG method, each edge is assembled once whereas for the FEQF method assembly is done element by element. Thus for the FEQF scheme, each edge is effectively assembled twice. The results from the advanced and simple model sets followed same assembly time trend. PAGE 82 82 An interesting difference between methods occurred when the mesh element nodes were displaced as a test of non ideal mesh conditions. With the exception of the g ate oxide channel interface and outside boundaries, each node inside the nMOS was randomly displaced by up to 40% of the initial grid spacing ( Figure 3 6 ) The randomization of the grid created a large number of negative edge coupling s which implies non Delaunay mesh elements throughout the structure. The negative coupling values were not zeroed out. For 2 D, the non Delaunay conditions were created by randomizing quad diagonal nodes. Non Delaunay conditions in 3 D simulations were cre ated using tetrahedron element types. Using the default N MOS structure (normal ideal grid) as a baseline, the results for both FEQF and FVSG methods were compared against equivalent structures with randomly displaced nodes. For both 2 D and 3 D simulations the FEQF method performed very accurately in terms of I ON for both normal and randomized grids ( Figure 3 7 ). However, the I ON results for the FVSG method deviated by a large amount, especially at small grid spacings. As the node ran domization was reduced, the FVSG method increased in accuracy. When using the FVSG method, solution convergence was a problem for the 3 D nMOS device simulations if non Delaunay elements were predominant. For both tetrahedra and bricks, if the mesh element s under the MOS gate were too flat (> 5:1 width:depth ratio) the FVSG solution would not converge. The FEQF method did not have trouble converging with this ratio. 3.4.2 Charge Collection simulations To examine the impact of the discretization methods on single event simulations, a reversed biased N+/P diode was used since it is a good representation of the source/drain junctions that are responsible for charge collection in MOSFETs. A charge PAGE 83 83 collection simulation was performed in 3 D since in 2 D, all quantities are assumed to be extruded into the third dimension which leads to a misrepresentation of the charge density. A 3 D N+/P diode was created as a template and tested with both tetrahedra and brick elements. A charge cloud based on the SPA equations in c hapt er 2 was generated into the depth of the device to model the electron hole pairs that are generated during a particle strike. Both methods converged well for DC bias conditions. However, the simulation results show that the FEQF method converged more reli ably in the transient domain than the FVSG method for different mesh spacings and charge concentrations. This could imply that the FEQF scheme handles isotropic current flow with more stability. This explanation is substantiated by the numerical stability problems that have been observed in the past for 3 D FVSG charge collection simulations [ 4 ]. In terms of mesh element types, both the FVSG and FEQF methods converged better for bricks than for tetrahedra, especially at high charge concentrations. A qualita tive comparison between the two methods is given in Figure 3 8 for single event upset simulations. Another major difference between discretization methods was noticed in their transient simulation times. The FVSG required more Ne wton steps for every solution time step than the FEQF method. Even though the assembly time for the FEQF method is ~22% longer, the total simulation time, on average, for a charge collection transient was less th an that of the FVSG method ( Figure 3 9 ). Because detailed charge collection simulations in 3 D often take a day or more to complete, this time savings could be quite significant. PAGE 84 84 3.5 Discretization Method Summary For the short channel nMOS simulations, both the FVSG and FEQF meth ods gave agreeing I ON results over a variety of grid spacings and element types. However, for a MOSFET grid with a non ideal mesh (non Delaunay elements), the FVSG method is prone to inaccuracy suggesting a high sensitivity to mesh alignment at channel int erfaces. Based on these results, the FEQF approach would most likely provide more accurate results for rough or curved interfaces or situations where meshing is non ideal. However, the FEQF method has the disadvantage of a longer DC simulation time due to a longer assembly time. For 3 D charge collection simulations, the FEQF method outperformed the FVSG approach due to a higher convergence rate which may be due to a better handl ing of isotropic current flow. The total transient simulation time was also le ss for the FEQF method. It should be noted that for 2 D charge collection simulations, both methods worked well and gave the same results with the exception that the FEQF gives a 1 3% reduced current transient peak. Even though the FVSG method is by far th e most accepted discretization scheme in practice today, the simulation results show that the finite element quasi Fermi discretization approach is a viable and in some cases preferable alternative for 3 D single event simulations PAGE 85 85 Table 3 1 Physical mo del sets used for the simulation comparisons Model Set Value s Simple Constant Mobility ( n,p = 150 cm 2 Advanced Philips Mobility, Lombardi Mobility, Velocity Saturation, SRH & Auger Recombination PAGE 86 86 Figure 3 1 Two dimension elements. A) Triangular. B) Rectangular. C) Pentagonal. Note that the rectangular and pentagonal elements have triangular equivalents. Figure 3 2 Three dimension elements. A) Tetrahedron. B) Hexahedron. C) Prism. Note that the hexahedra elements can be divided into tetrahedral/prism equivalents. Figure 3 3 Two dimension example for an area A i associated with a node (represente d by circles) for generalized box discretization. PAGE 87 87 Figure 3 4 NMOS I ON currents using the advanced physical model set for a variety of grid spacings. Figure 3 5 Percent change in solution time per Newton step for the FEQF when compared to the FVSG (orange baseline). Based on the NMOS template with quad diagonal elements and advanced physical models. PAGE 88 88 Figure 3 6 An example of a perturbed mesh for the NMOS simulations. Note that current flow is not aligned with the grid in the channel region. Figure 3 7 The FVSG method loses accuracy for highly non Delaunay mesh conditions in the NMOS channel. The FEQF method is less affected by the non ideal mesh conditions. Average based on 10 simulations per grid spacing. PAGE 89 89 Figure 3 8 Qualitative s ingle event upset comparis on for both discretization methods. PAGE 90 90 Figure 3 9 Normalized average total transient simulation time. The average was taken over 15 simulations each with difference charge concentrations. PAGE 91 91 CHAPTER 4 DEVICE GRID AND BOUNDARY SCHEMES 4.1 Introduction This chapter builds on the discussion of discretization methods in the previous chapter by closely examining the grid generation and boundary edges for simulation devices used for single event simulations. The focus of this ch apter is on finding ways to reduce simulation time, since SEE simulations are very time intensive. This chapter discusses two new proposed methods that offer simulation time savings while maintaining a high level of accuracy in results. The first section w ill describe an adaptive gridding scheme which reduces the number grid points (and simulation time) in real time for a single event transient. The second section will discuss a proposed diffusive boundary scheme that can be used for the non contacted outer edges of a simulation structure. 4.2 Adaptive Gridding Continued advancements in technology computer aided design (TCAD) and physical modeling have enabled increasingly complex device structures to be characterized. However, radiation effects simulations int roduce additional complexities that modern TCAD tools are not well designed for. For SEU, the grid generation (a.k.a. mesh generation) of the device structure is a key area in need of improvement. SEU simulations introduce great complexity since a high gri d density around the strike region is required to resolve the carrier movement from the electron hole (e h) pairs generated by the particle strike. This requires the SEU modeler to create a customized grid for each device simulation structure and account f or variables such as particle strike path, energy, and angle of incidence. However, once the simulation has started and the PAGE 92 92 transient progresses, the particle strike induced carriers diffuse widely throughout the device and a dense grid in the strike regio n is no longer needed. Because the total simulation time is directly proportional to the number of grid points, adapting the grid to the needs of the simulation in real time during SEU transients could result in significant time savings. Although grid ada ptation techniques suitable for steady state simulations have been developed, the SEU modeler needs the capability to dynamically adapt the grid as the transient progresses, as discussed by Dodd in [ 4 ]. To our knowledge, no such transient gridding techniqu es exist in currently available TCAD tools. In addition to SEU modeling, adaptive gridding would be beneficial for the general purpose simulation and characterization of modern devices. For instance, the 2007 International Technology Roadmap for Semiconduc tors (ITRS) states that advances in grid adaptation are a priority for TCAD tool development since devices are becoming increasingly complex [ 39 ]. This work proposes a practical way to adaptively refine and coarsen the grid around a strike path during a S EU simulation [ 40 ]. The adaptive grid scheme reduces the time spent by the SEU modeler on customizing the grid and reduces the total simulation time while maintaining a high level of accuracy in results. As with the rest of this work, the TCAD tool used fo r this work is FLOODS since the code is readily customizable [ 1 5 ]. The adaptive grid algorithms in this work were implemented in FLOODS using Tcl/Tk and C++. 4.2.1 Minimizing Discretization Error In order to have a basis for grid refinement and coarsening, the k ey mechanism for minimizing discretization error for single event simulation need to the examined. For PAGE 93 93 SEE, the key parameter is the charge generation and collection. Therefore, the discretization error relating to the charge should be minimized in order f or the simulation results to be accurate. Building on the discussion in the previous chapter, consider ( 1 33 ) The volume integral is approximated by using the value associated with node i and multiplying it by the area as shown in Figure 3 3 The error expression for this approximation is straightforward and can be writ ten as ( 1 37 ) where the error is proportional to the grid spacing squ ared and the second derivative of the charge the second derivative of potential is equal to the charge. Thus, to minimize discretization error f or SEU simulations, the grid spacing should be very small in depletion regions and wh ere the charge is high (i.e. strike path). 4.2.2 Simulation Time Tradeoff The grid spacing is proportional to the discretization error as discussed in [ 31 ] More specifically, a smaller distance between grid points (nodes) results in a smaller discretization er ror and thus a more accurate simulation result. However, an interesting tradeoff between accuracy and solution time exists. The solution time increases rapidly with the number of grid points as ( 1 38 ) where the term varies between 1.5 and 1.75, t is the solution time and m is the number of grid points [ 32 ] Figure 4 1 illustrates the dependence of solution time on the number of grid points. In this example, a two dimensional uniformly doped resistor with PAGE 94 94 dimensions of 1.01.0 m was simulated and the number of grid points within the device was varied. As the number of grid points increases, the solution time quickly increases as well. This illustrates the need to carefully limit the number of grid points so that the simulation time does not become excessive. Unfortunately, SEU simulations are ty pically very time intensive for two reasons. First, a high number of transient solution steps are required to simulate the current voltage response and the charge collection process. Second, in addition to the gridding required for the steady state solutio n, more grid points are needed around the region of the particle strike so that the numerical solution converges with accurate results. 4.2.3 Adaptive Grid Scheme Methodology The proposed adaptive grid scheme is based on creating a set of individual grids with v arying levels of complexity (grid points). A flowchart of the scheme is given in Figure 4 2 First, a grid is generated that is suitable for running standard steady state DC simulations. A DC mesh should be refined around the depletion regions, junctions, and interfaces (i.e. MOSFET channel) to reduce discretization error. Next, the steady state bias conditions are simulated for the device and then the electron hole pair distribution is generated to model the particle strike. For this work, the e h pair pro files are given as a constant for the transient simulation and are based on the models described in the next section. Once the e h pair distribution is known, the grid is refined a user specified number of times around the strike region. To achieve accurat e SEU simulation results, a device structure with a very coarse initial grid in the bulk region will require more refinements than a structure that already has a very refined grid. The refinements are based on evaluating the boundaries of C ref where PAGE 95 95 ( 1 39 ) where n and p are the electron and hole densities. Because C ref is a good approximation of the electron hole pair distribution, it allows for a very straightforward refinement of the strike path, as shown later in Figure 4 8 and Figure 4 10 If refinement were based only on the electron or hole density, any heavily doped region (i.e. MOSFET gate, source/drain) would be refined further and possibly unnecessarily so. Although the net charge of an electron hole pair is zero, the resulting separation of carriers (e.g. funneling, drift, diffusion) determines the charge collect ion. If the area around the strike path C ref is poorly gridding, a large amount of discretization error is introduced as in (6), where the charge discretization error is a function of grid spacing. Thus, gridding around C ref insures that any discretization error in approximating the e h charge distribution is minimized. The grid refinement process works by taking a specified region of the grid and then dividing each grid element wi thin the region. For instance, rectangular grid (quad) and triangular element s will be split into four smaller elements Figure 4 3 illustrates the refinement of a 11 m structure made up of quad elements. In this example, three refinements are performed on a Gaussian function that is similar to an ion strike track. The Gaussian h as a 1/e radius of 50 nm, a peak e h concentration of 1.110 20 cm 3 and the grid is refined inside the 10 18 10 19 and 10 20 cm 3 contours of C ref After every refinement, each grid is stored so that a collection of different grids is accessible to the devi ce simulator for later use ( Figure 4 2 ). In FLOODS, the grid generation/storage algorithm is fully adaptable to rectangular and non rectangular elements (i.e. Delaunay triangular mesh). Thus, there is no gain from using one element PAGE 96 96 type over another from t he standpoint of grid generation/storage efficiency. The simulation structures in the next section are built using rectangular elements where the five point elements are divided into triangular elements. Since the grid should be aligned in the direction of current flow under steady state conditions to minimize discretization error [ 36 ], rectangular elements work well since current flow in MOSFETs and diodes is laminar in nature. Following the adaptive refinement around the strike path Cref, the transient simulation is started as in Figure 4 2 As the transient simulation progresses, the grid is continually coarsened and eventually resolves back to the original grid used for the DC solution. For each refinement or coarsening step, the values for every simul ation variable (i.e. electron and hole density, doping, electrostatic potential) are interpolated from one grid onto another. In this work, the grid is coarsened each time the peak carrier density in the strike region falls by an order of magnitude. This e nsures that the grid coarsening process does not occur until the charge has started diffusing throughout the device. If the grid is coarsened too soon, valuable information on the charge distribution is lost. When the grid is coarsened, it is inevitable t hat some error gets introduced when interpolating the variables (i.e. potential, carrier density) from one grid onto another. To compensate for the new grid and associated variables, the simulation tool needs to return to a small time step in order to damp en the error that was just introduced. FLOODS self estimates each time step and uses the TR BDF time discretization method [ 41 ]. For the adaptive grid simulations in the next section, the first time step after coarsening typically fell into the femto and picosecond range. Although the grid PAGE 97 97 coarsening process adds time steps, the benefit from having less grid points still results in an overall time savings. However, in the case of a grid being coarsened too many times, the resulting addition of time steps w ould start to negate the benefits using of a coarse grid. In this work, the grid is coarsened a maximum of three times during the transient simulation. 4.2.4 Simulation Results Single event transient (SET) simulations were performed to compare the results obtain ed using the adaptive grid scheme versus two different static grid schemes. The first static grid scheme uses a uniform grid over the entire structure. An ultra dense uniform grid will yield the best simulation result, but the longest simulation time. The second static grid scheme uses a grid that has been refined around the junctions and the particle strike region, similar to a customized grid that an experienced SEU modeler might create. It is important to note that the customized grid requires some TCAD experience and may not be an option for an inexperienced TCAD user in SEU modeling. Two different sets of simulations were run to compare the grid schemes. For the first simulation set, the grid schemes are compared for a laser induced current transient in a reverse biased N+/P diode structure, similar to the scenario described in [ 42 ]. For the second simulation set, a particle strike path is generated in an nMOSFET device with a 90 nm gate. Each simulation uses the Philips unified mobility model, the Shock ley Read Hall and Auger band to band recombination models [ 37 ]. An overview of the simulation variables is given in Table 4 1 It is important to note that a comparison to experimental data is neglected since the focus this work is to examine the tradeoff between solution time and discretization error from the grid. Furthermore, since an PAGE 98 98 increasingly high grid density converges towards a specific result, the results from the structure with the most grid points can be viewed as the best possible result (smal lest discretization error). FLOODS is currently limited to adaptive gridding in 2 D but the results would also be applicable for 3 D applications. Arguably, the easiest particle strike to grid would be Gaussian in form, uniform in depth, and normally incid ent to the surface. However, to illustrate the benefit of the adaptive grid scheme, each simulation set uses a unique e h pair profile that is more challenging to grid. For the N+/P diode simulation set, the number and distribution of e h pairs generated by a laser pulse is calculated by using the single photon absorption equation developed by McMorrow [ 22 ]. This model is discussed in chapter 2 and is given by equation ( 1 7 ) For the second simulation set, the generated electron hole pair profile is modeled using an angled cylindrically symmetrical Gaussian profile. The Gaussian profile had a 1/ e radius of 5 nm, terminated at a depth of 0.4 m, a LET value of 0 .1 MeV cm 2 /mg and an incident angle of 30 degrees. Figure 4 4 shows the carrier dist ribution for the SPA model and the Gaussian profile. The peak carrier concentrations for the SPA and Gaussian profiles were 8.8410 18 cm 3 and 8.2110 19 cm 3 respectively. 4.2.4.4 N+/P Diode Simulation For the first simulation set, single event transient simulations for an N+/P diode were run to compare the results obtained using the adaptive grid scheme versus a customized and uniform grid scheme. The diode simulation structure is 30 30 m and consists of a heavily doped n+ region (10 20 cm 3 ) in a p well (10 18 cm 3 ) that resolves into a p type substrate (10 16 cm 3 ). The n+ and p PAGE 99 99 electron hole pairs for the diode is shown in Figure 4 4 The simulation results for the current and collected ch arge versus time are given in Figure 4 5 and Figure 4 6 As the uniform grid was coarsened, the e h charge profile interpolation error was increased and the charge was overestimated. Additionally, with uniform grid coarsening, the depletion region was over estimated which resulted in a higher collected charge ( Figure 4 5 ). As expected, the customized grids had a shorter simulation time than the uniform grids, shown in Figure 4 7 However, only the ~8,000 and ~15,000 point customized grids had the same accura cy as the ultra dense uniform grid. The adaptive grid scheme was simulated for different numbers of refinement levels and X =3 as found to give the best results in terms of time savings. The simulation time versus the collected charge is given in Figure 4 7 and illustrates the importance of coarsening the grid in real time as the SET progresses. The adaptive grid scheme strikes a good balance between the simulation time and accuracy in collected charge. For example, a diode structure with a uniform grid of ~23,000 points would take more 10 times longer to simulate than the adaptive grid for the same result. Comparing Figure 4 7 and Figure 4 4 it can be seen that the areas of highest grid refinement correspond to areas of highest e h pair density. For both t he diode and N MOS simulations, it was found that refinement worked the best when starting around the C ref contour of 10 15 cm 3 Refinement at this C ref value covers the outer boundary of the strike path and limits discretization error from potential contou r deformation and diffusion as the transient progresses. For example, for the adaptive grid at X =3, the PAGE 100 100 refinements were done about the C ref contours of 10 15 10 16 and 10 17 cm 3 As a side note, the simulations were performed using a 2.93GHz quad core pro cessor and a normalized time of 100 and 1000 in Figure 4 7 represents a simulation time of 70 and 700 minutes respectively. 4.2.4.5 N MOS Simulation For the second set, SEU simulations for an nMOSFET were performed to compare the results obtained using the adaptive grid scheme versus a customized and uniform grid scheme. The nMOSFET simulation structure is based on the 90 nm technology node with a bias of 1 V applied to the drain. The oxide thickness was 2 nm and the physical gate length and height were 90 nm and 60 nm respectively. The doping profile was based on analytic functions and values given in [ 43 ] [ 44 ]. The nMOSFET simulation time versus the number nodes is given in Figure 4 9 and further illustrates the benefit of adapting the grid in real time during the transient. The uniform grid with ~52,000 grid points was used as the baseline for the collected charge since it had the highest grid density. As the uniform grid was coarsened, the interpolation error increased and the e h charge profile was overestimated as with the N+/P diode case. Likewise, the depletion region was overestimated due to low grid densities and resulted in a higher error in collected charge. The customized grids had a shorter simulation time than the uniform grids and the 15,000 point cust omized grid had the same accuracy as the ultra dense uniform grid. Again, the adaptive grid scheme finds a good balance between the simulation time and accuracy in collected charge for X =3. For example, the nMOS structure with a uniform grid of ~23,000 ta kes about 10 times longer to simulate than the adaptive grid for the same result. The adaptive grid for the NMOS simulations is shown in Fig. 4 9 at PAGE 101 101 X =3 levels of refinement. Comparing Figure 4 10 and Figure 4 4 it can be seen that the areas of highest gr id refinement correspond to areas of highest e h pair density. In this example, the strike path did not cross any insulator boundaries (i.e. STI, gate oxide). However, a strike path that traverses through insulators should be refined to minimize the discre tization error of the electrostatic potential. 4.2.5 Adaptive Grid Summary This section presented an adaptive grid scheme for SEU simulations with results that show the proposed scheme can offer significant simulation time savings while preserving accuracy. The time saving benefits of the proposed scheme would be especially useful for the automation of SEU simulations. Programs such as Monte Carlo radiative energy deposition (MRED) are used to generate very large numbers of individual single event descriptions fo r 3 D structures such as latches and SRAM cells [ 45 ]. A program like MRED could generate the electron hole pair charge distribution and then use a device simulation tool with the adaptive gridding scheme to simulate the each SEU automatically [ 46 ]. This wo uld eliminate the need for an experienced TCAD user to have to custom grid each simulation structure. An additional benefit of the proposed scheme is that the refinement parameters can be adjusted by the user to yield more accurate results (denser adaptiv e grid) or a faster simulation time (coarser adaptive grid). Although the results are only shown for 2 D simulations, the adaptive grid scheme could be applied to 3 D simulation structures where the time benefit may be even greater since larger differences in grid density occur. PAGE 102 102 4.3 Boundary Sinks The outer edges of a device simulation structure that are not associated with contacts (i.e. ohmic, schottky) are important for single event simulations. First, the outer boundary placement affects the simulation time A larger device boundary usually contains more grid points, which in turn, increases simulation time. Second, the device boundary affects the accuracy of the simulation results. If an outer boundary is too small, it will affect the key operating regions of the device and will lead to inaccurate results. Normally, one can define a reasonably small boundary for a device. For example, an NMOS device will not need to have the entire source/drain or substrate under the p well defined in order for the results t o be accurate since the channel determines the current output. However, with single event effects, carriers diffuse widely throughout the device requiring a larger outer boundary to be created. If the boundary is too small, charge collection will be overes timated since most TCAD tools use reflective boundaries (zero flux condition) at the non contacted device edges. The work in this number of electrons and holes to cross a non c ontacted boundary. This allows for the approximation of a larger device outer boundary than what will actually be created (thus a simulation time savings). In the following subsections, the proposed boundary sink with respect to different device boundary s izes will be discussed. 4.3.1 Boundary Condition Overview For most semiconductor device problems, both Neumann and Dirichlet boundary conditions occur for the PDEs (i.e. Poisson and continuity equations) [ 36 ]. In FLOODS and most device simulation tools, a non c ontacted boundary is a boundary in which no flux is allowed to pass. This condition for the flux F (i.e. E, J n J p ) can be written as PAGE 103 103 ( 1 40 ) where n is the unit normal vector to the contour of integration as discussed in the previous discretization chapter. This condition is referred to as a homogenous Neumann boundary c ondition. This boundary is simple to implement since it means that the integration along the boundary edges is completely ignored. However, the homogenous Neumann boundary induces reflective symmetry. For example, consider a diode with only the left half o f the device versus the full device. Figure 4 11 shows the potential contours of the reverse biased diode where it can be seen that when the diode is cut in half, the solution for both devices is the same. In other words, if the current from the half diode were multiplied by two, it would give in the same result as the full diode. The carriers are reflected at the boundaries. However, if the outer device boundaries are large enough, the carriers will recombine before any reflective boundary issues impact the results (at the cost of more grid points and longer simulation time). Dirichlet boundaries are edges for which the solution variables ( n p ) are fixed for the PD Es and are typically used for contacts. For example, for an n type ohmic contact, the electrostatic potential is fixed at the boundary as ( 1 41 ) where V applied is the applied potential at the contact. Furthermore, in equilibrium, the quasi he applied potential as ( 1 42 ) PAGE 104 104 For an ideal ohmic boundary condition or contact, there is no limit on the amount of current flowing through the contact interface (a.k.a. infinite surface recombination velocity). To define the new proposed diffusive boundary sink, a few modifications are made to the Neumann and Dirichlet co nditions. First, the sink is formulated so that the This leads to the homogenous Neumann condition of ( 1 43 ) where E is electric field. In other words, the potential has reflective symmetry at the boundary sink edge. This conditi on is employed for the potential since the proposed boundary sink should not behave like contact of any form or function (i.e. supply, ground), as it would affect normal device operation. Next, consider the electron flux at the boundary sink (the following arguments apply the same way for holes). For a homogenous Neumann condition, the electron flux would be zero. However, to approximate a larger boundary or device volume, some ple, Fi gure 4 12 whi ch compares both boundary types. In equilibrium, the electrons should not be freely flowing past the boundary sink. However, for a particle strike, the high concentration of electrons diffusing throughout the device generates an excess of elec trons at the boundary sink with respect to equilibrium levels. One way to allow electrons past the boundary is to assume a finite surface recombination velocity that is a function of the diffusion length. This can be formulated as ( 1 44 ) PAGE 105 105 ( 1 45 ) with ( 1 46 ) and ( 1 47 ) where n and p are the electron and hole densities, n eq and p eq the equilibrium densities U s n and U s p the surface recombination rates, D n,p the diffusion coefficient, L n,p the diffusion length, n,p the carrier mobility, n,p the carrier lifetime and x the distance from the boundary. For equations ( 1 44 ) and ( 1 45 ) the diffusion terms work out into units of a finite recombination velocity. As a side note, the carrier lifetime is a function of doping and temperature as discussed in chapter 2. Also, the carrier mobility terms can be written as functions of doping leve ls, carrier concentration and temperature as will be discussed in chapter 6. Therefore, the recombination rates for the proposed boundary sink approach are functions of doping, carrier concentration, and temperature. The next section will show single event simulation results using the both the reflective and diffusive boundary sink conditions. 4.3.2 Simulation Results To test the diffusive boundary sink, a simple reverse biased N+/P diode was used (similar to Chapter 2). The initial 2 D simulation structure was 3 0 40 m in width and depth as in Figure 4 13. To mimic an ion strike, the electron hole distribution was modeled using equation ( 1 3 ) and correlates to a c onstant LET of 10 MeV cm 2 /mg. The peak carrier concentration of the strike is 8.2110 1 9 cm 3 has a 1/ e radius of 50 nm and PAGE 106 106 terminates at a depth of 30 m. The physical models used were the UF high injection mobility model (Chapter 6) and the Auger and SRH recombination models. For the simulation comparison, the device width was varied to the values 10, 30, 100 and 200 m. For the first simulation set, standard reflective boundaries were used for the right and left edges of the device. The results for each width are shown in Figure 4 14 and the trend is that as the width is decreased, an increase in collected charge occurs. Even with recombination, for the reflective boundaries for small widths, an excess of charge collection is observed. However, for exces sively large widths (100 and 200 m) the results converge to a specific answer since most of the particle strike induced carriers have recombined by the time they reach the boundary. Note that the charge collection deviations in Figure 4 14 correspond to the diffusion component of charge collection (t > 10 8 s). Very little change in drift/funneling current collection (t < 10 8 s) is observed since drift/funneling is more of a function of the depletion/funnel field region than the outer boundaries. Thus, c urrent transient plots are not shown for this section. When using diffusive boundary sinks, the error in total collected charge is well controlled. As shown in Figure 4 15, results for the device with boundary sinks converge to the same result as the larg e device with reflective boundaries at 200 m. Interestingly, the 10 m wide device with boundary sinks converges to the same total collected charge as the reflective 200 m wide device. Note that for the 10 m device width, the smaller boundary causes mor e charge to be channeled toward the top contact (e.g. more electrons diffuse into the depletion region) around the time of 5 10 8 seconds in Figure 4 15. However by formulation, the diffusive boundary sinks (and PAGE 107 107 surface recombination rates U s ) are a funct ion of excess carriers (i.e. n n eq ) which compensates for this channeling effect at the smaller widths. A direct comparison of both boundary types is given in Figure 4 16. When the reflective boundary device width is reduced, an error in collected charge occurs. However, the device with diffusive boundary sinks is immune to this effect. Figure 4 17 adds a simulation time comparison, where is clear that there is a significant time savings benefit to using the boundary sinks. 4.3.3 Boundary Sink Summary A propose d diffusive boundary sink approach was formulated and shown to give excellent time savings results. Although it allows the TCAD user to reduce the outer boundary size for SEE simulations, care should still be taken in choosing the device boundary edges. Fo r example, the boundaries should never be reduced to the point of affecting the steady state operation of the simulation device. Additionally, the outer edges of the simulation structure should always surround the particle strike path. PAGE 108 108 Table 4 1 Overvie w of the adaptive grid simulation variables Simulation Set Set 1 Set 2 Structure Type N+/P diode N type MOSFET Generated electron hole pair profile Single Photon Absorption [12] Energy = 13.5 pJ Cylindrical Gaussian LET = 0.1 MeV cm2/mg PAGE 109 109 Figure 4 1 An increase in m grid points results in an increase in solution time. Figure 4 2 Flowchart of t he proposed adaptive grid scheme. PAGE 110 110 Figure 4 3 Example of grid refinement on a Gaussian function. A B Figure 4 4 Electron hole pair distributions used in simulations. A) Single ph oton absorption, laser energy = 13.5 pJ, radius = 2 m. B) Cylindrical Gaussian, PAGE 111 111 Figure 4 5 N+/P diode 2 D simulation results comparing current transients for the uniform and a daptive grid schemes. Figure 4 6 N+/P diode 2 D simulation results comparing collected charge versus time for the uniform and adaptive grid schemes. PAGE 112 112 Figure 4 7 N+/P diode results. The number of grid points is given next to each data point. The results were normalized and a value of 1 on both scales represents the lowest discretization error (y axis) and the fastest simulation time (x axis). Figure 4 8 Adaptive grid at peak refinement ( X =3) about the C ref contours of 10 15 10 16 and 10 17 cm 3 for the N+/P diode simulations. Grid points m=7,854. PAGE 113 113 Figure 4 9 nMOSFET results. The number of grid points is given next to each data point. The results were normalized and a value of 1.0 on both scales represents the lowest discretization error (y axis) and the fastest simulation time (x axis). Figure 4 10 Adaptive grid at peak refinement ( X =3) about the C ref contours of 10 15 10 16 and 10 17 cm 3 for the nMOS simulations. Grid points m=8,114. PAGE 114 114 A B Figure 4 11 Example of reflective symmetry using FLOODS A) Half diode cross section. B) Full diode cross section. Figure 4 12 Illustration comparing: A) homogenous Neumann boundary. B) proposed diffusive boundary sink. PAGE 115 115 Figure 4 13 Illustration showing the simulation structure with diffusive boundary sinks for two different widths. Figure 4 14 Collected charge versus time for a reversed biased N+/P diode with reflective boundaries. PAGE 116 116 Figure 4 15 Collected charge versus time for a reversed biased N+/P diode with diffusive boundary sinks Figure 4 16 Comparison of boundary types with respect to device width and collected charge. PAGE 117 117 Figure 4 17 Comparison of boundary types with respect to device width collected charge and simu lation time PAGE 118 118 CHAPTER 5 IMPACT OF STRAINED SILICON ON SINGLE EVENT EFFECTS 5.1 Motivation and Background develop innovative new processing techniques. Recently, a large amount of focus has been on using front end process induced stress to improve channel mobility and thus transistor I ON performance [ 47 ]. For the 45 nm technology node, the channel stress is induced using SiGe source/drain implants and compressive capping layers for PMOS devic es and tensile capping layers for NMOS devices [ 48 ]. In order to accurately characterize single event effects for modern CMOS transistors, the impacts of strained silicon technology need to be considered. Although CMOS devices with feature sizes 22 nm and smaller have been reported, the 45 nm node will be the focus of this chapter. The reason is that newer device technologies are not implemented in spaceborne systems until they have been well characterized in terms of both single event and total ionizing d ose response. As of the year 2010, the 45 nm and 65 nm nodes have been the focus of much experimental SEE work in the radiation effects community. In this chapter, a brief overview of strained Si physics is given with respect to electron and hole mobility Next, an overview of stress and strain tensor matrixes is given in order to better understand piezoresistance. Following the discussion of stress, a piezoresistive mobility model that is function of crystallographic orientation is formulated. Although, t his model is currently used in other modern TCAD simulation tools, the crystallographic dependencies are not described or universally formulated for these tools [ 21 ]. Following the piezoresistance modeling discussion, the piezoresistive mobility model is c ompared against experimental results for a uniaxially strained N+/P PAGE 119 119 diode, where it is shown that the results match well. Then, single event transient predictions are made for strained silicon CMOS devices at the 45 nm node. For these simulations, FLOOPS i s used to calculate the process induced channel stress in all directions (i.e. longitudinal, transverse). In the previous chapters, the focus was on improving the simulation tool in terms of discretization, gridding and boundary methods. However, for this and the remaining chapters, the main focus will be on the physical modeling of carrier mobility in silicon. 5.2 Brief Overview of the Physics of Strained Silicon Although the effects of strained silicon on mobility have been studied for many years, it has been the topic much interest for the past last decade since can be used to enhance MOSFET channel mobility. Mobility in silicon it is commonly expressed in a generalized form as ( 1 48 ) where m is the mean free time between collisions (1/ is the scattering rate) and m* the conductivity effective mass [ 49 ]. The effective mass and scattering terms in silicon are changed by stress. The remainder of this section gives brief overview of the physics behind strained silicon for both electrons and holes. A much more thorough overview of the physics of strained Si is given by Sun et. al [ 50 ]. For the case of electrons, strain induced mobility enhancement is best explained by describing the conductivity effective mass and scattering. Figure 5 1 shows the conduction band for bulk unstrained Si at room temperature. The conduction band is 6 ) where the degeneracy reflects the cubic symmetry of the Si lattice [ 51 ]. However, the effective mass of each ellipsoid is PAGE 120 120 anisotropic and the longitudinal mass m l (parallel to axis) is larger than t he transverse mass m t (perpendicular to axis). The electron conductivity mass m* for unstressed bulk Si can be written as ( 1 49 ) where m 0 is the free electron mass, m l =0.98 m 0 and m t =0.19 m 0 [ 51 ] For the case of a device on a (001) wafer, advantageous strain sp 6 4 (in plane) and 2 (out of plane) valleys as in Figure 5 1. 2 valleys means that they are preferentially populated by electrons and the electron mobility improves due to a reduced in plane effective mass m* Additionally, it is believed that intervalley phonon scattering is reduced due to the splitting of the conduction valleys [ 52 ]. Revisiting equation ( 1 48 ) it can be seen that strain can be used reduce scattering and conductivity mass for electrons which results in an increase in mobility. For the case of a NMOS channel, the 2 4 valleys are already split due to the gate bias. Thus, the main contribution to electron mobility enhancement is likely due to scattering (i.e. phonon, surface roughness). For holes, an increase in mobility relates to the valence band warping. T he valence band structure is more complex than the conduction band for both unstrained and strained Si. For unstrained Si at room temperature, holes occupy the heavy and light hole bands at the top of the valence band. When strain is applied, the hole effe ctive mass becomes highly anisotropic due to band warping. Subsequently, the valence energy levels breakup into separate heavy, light, split off bands [ 51 ]. Analogous to electrons, holes preferentially occupy the top band at higher strain due to the strain PAGE 121 121 induced energy level splitting as in Figure 5 2 and experience a lower in plane mass. A high density of states is required to sufficiently populate the top band and it has been found that uniaxially compression in the channel direction <110> for (100) and (110) wafers gives desirable results [ 51 ]. Additionally, as stress levels greater than 1 GPa are induced, hole intervalley scattering is reduced, resulting an increase in hole mobility. For modern CMOS transistors, uniaxial stress along the channel direct ion <110> is used to enhance mobility for both NMOS (tensile stress) and PMOS (compressive stress) devices. 5.3 Piezoresistance Mobility Model Many models exists for estimation the change in mobility due to stress. Of these models, the piezoresistance mobility model is the most computationally efficient (important for SEE) and practical for device simulations. In this work, FLOOPS is used to calculate the front end process induced stress profiles for the single event N+/P diode and MOSFET simulations described later in this chapter. An understanding of how stress is calculated is necessary to describe how stress is used as an input to the piezoresistance mobility model. This section starts with a discussion on linear elasticity and then describes how the strain and stress tensor matrixes are formed. 5.3.1 Linear Elasticity Linear elasticity is a property of solid materials that determines how objects deform and become internally stressed due to externally applied loading conditions. The near elasticity is that a linear relationship between strain and stress exists between the corresponding axis components of stress and strain for conditions that do not produce yielding (permanent deformation). This assumption is PAGE 122 122 commonly used for finite e lement analysis of structures such as semiconductor devices [ 53 ]. is directly proportional to the external load (as long as the load does not surpass the elastic limit ). In one ( 1 50 ) where F is the restoring force exerted by the material, k is the stiffness associated with the material, and x is the displacement of the end of the material from its equilibrium position [ 54 ]. The stiffness k is a measure of how resistant the material is to ex ternal forces. In a process simulation tool such as FLOOPS, a stiffness matrix k is used to example, a 1 D spring element associated with two nodes is shown in Figure 5 3. The relationship between nodal forces and nodal displacements shown in Figure 5 3 can be written in matrix form as ( 1 51 ) where k ij are the element stiffness coefficients of the k matrix and d the associated nodal unit displacements in the x axis. If a force is applied to the spring, an equal and opposite force is generate d. x (or strain) related to equation ( 1 50 ) For instance, if the nodes in Figure 5 3 are subjected to tensile forces, the spring will deform by expanding and the d 1 d 2 displacement values in equation ( 1 51 ) will change. At a higher level (i.e. device grid), a large number of elements (and nodes) exist and a global stiffness matrix needs to be assembled such that PAGE 123 123 ( 1 52 ) This work uses the common assumption that silicon acts as a linear elastic material for the stress inducing processing conditions that are common in modern CMOS technologies. 5.3.2 The Strain and Stress Tensors Strain ( ) is a unitless param eter that relates to the deformation of a solid body that is subjected to a force. It is equal to the change in length in a given direction divided by the initial length L simply as ( 1 53 ) or in terms of the normal components ( 1 54 ) where u v, and w represent the displacements in the x y and z directions respectively. For linear elastic materials like silicon, the cross section becomes narro wer when strain and is written as ( 1 55 ) in the x direction with respect to a change in y plus the displacement in the y direction with respect to a change in x [ 53 ]. For exa mple, the shear strain xy can be written as PAGE 124 124 ( 1 56 ) All nine norm al and shear strain components can be combined in a strain tensor matrix ij as ( 1 57 ) A tensor is an object which extends the idea of scalar, vector, or matrix, and does not vary from the transformations of coordinates. In static equilibrium, the shear component are equal (i.e. xy = yx ) and the strain tensor can be condensed into six components as ( 1 58 ) Stress follows a similar fo rm to that of strain in terms of matrix formulation. Stress S ) within a deformable body as ( 1 59 ) As with strain, stresses have normal and shear components. Shown in Figure 5 4 is a three dimensional infinitesimal element in Cartesian coordinates. Normal forces act perpendicular (normal) to the faces and shear forces ac t along each face of a body. Tensile forces are positive and compressive forces are negative. Similar to strain, stress PAGE 125 125 ij is ( 1 60 ) In static equilibrium, some of the shear stresses are equa l by symmetry (i.e. xy yx ) and the stress matrix can be reduced to ( 1 61 ) In a line arly elastic material, the stress is linearly proportional to the strain. Using the ( 1 62 ) where D is the stress/strain matrix (or constitutive matrix) [ 54 ]. Since silicon can be approximated with isotropic elastic properties, the constitutive matrix can be written as ( 1 63 ) PAGE 126 126 where E addition to a derivation of the D matrix, a detailed description on how strain and stress is coded in FLOOPS is given by Randall in [ 53 ]. 5.3.3 Piezoresistance Definition Now that the stress tensor matrix has been defined, the piezoresistance model can be desc ribed. Piezoresistivity is the change in electrical resistivity ( ) due to mechanical stress ( ) It involves the relationships, both linear and nonlinear, between the electric field E i current density J j and mechanical stress kl [ 55 ]. The change in ele ctric field dE i with stress and current can be expanded in a McLaurin series as ( 1 64 ) w conductivity. The first term ( dE i / dJ j ) is the electrical resistivity ij a second rank polar tensor. The fifth term ( d 2 E i )/( dJ j jk ) is the fourth rank polar tensor that describes the dependence of electrical resistivity on stress [ 55 ] The odd rank polar tensors disappear in the McLaurin series since silicon and germanium are from a centrosymmetric (m3m) point group resulting in ( 1 65 ) Integrating both sides gives ( 1 66 ) where the stress induced change in resistivity is PAGE 127 127 ( 1 67 ) For the m3m point group there are three independent tensor coefficients. These three independent piezoresistance coefficients and are given by equations ( 1 68 ) ( 1 69 ) and ( 1 70 ) The coefficients can be reduced as such b ecause the stress tensor kl is symmetric in silicon thus k and l can be interchanged. Likewise, i and j are interchangeable because the resist ivity ij and strain ij tensors are symmetric. However, i and j cannot be interchanged with k and l It is important to note that the re lationship between the matrix and tensor coefficients involves a factor of two whenever ij is defined by i = 1 6, j = 4 6. For example 66 1212 and 44 1313 as shown in the following equations as ( 1 68 ) ( 1 69 ) ( 1 70 ) All other tensor coefficients are zero for silicon and other point group m3m crystals. It is common convention that the ijkl coefficients, the i j kl pairs c an be replaced as the following ( 1 71 ) PAGE 128 128 For equation ( 1 71 ) it is helpful to visualize the orientation as in Figure 5 5, where the z axis [001] often corresponds with depth into the device and the [110] orientation is in the same direction as a CMOS channel. For the general case of a tryclinic crystal, the shortened matrix form of the piezoresistance matrix would be ( 1 72 ) However, since silicon, germanium, and other crystals with cubic symmetry have only three independent tensor coefficients the previous matrix reduces to ( 1 73 ) Since the CMOS channel is in the <110> direction, equation ( 1 73 ) needs to be transformed. The fully transformable piezocoefficent matrix is given by the following equation as ( 1 74 ) PAGE 129 129 The full deriv ation of this matrix is given the Appendix. From the derivation, it is shown that the transformable piezoresistance coefficients in equation ( 1 74 ) are given by the following equations ( 1 75 ) ( 1 76 ) ( 1 77 ) ( 1 78 ) ( 1 79 ) ( 1 80 ) ( 1 81 ) ( 1 82 ) ( 1 83 ) where the l, m, n values represent directional cosine transformations given by ( 1 84 ) For the coordinate system in Figure 5 5 is axis and for y axis. For the case of a CMOS device with a channel orientation of [110], a value of =45 degrees should be use d for equation ( 1 84 ) PAGE 130 130 5 1) and are commonly used to consider mobility enha ncement under mechanical stress in silicon [ 60 ]. In addition to stress, piezoresistance is also a function of impurity concentration and temperature as shown by Kanda [ 61 ]. The dependence of the piezoresistance on impurity concentration and temperature can be written as ( 1 85 ) where E F is the Fermi level, F 0 is the Fermi in tegral of the order 0, and F 0 the first derivative. The term P ( N,T ) is commonly referred to as the piezoresistance factor. An example of the piezoresistance factor for n type silicon is given in Figure 5 6 where it is evident that as temperature and impur ity increase, the piezoresistance effect decreases. It should be noted that although the piezoresistance factor is a function of impurity concentration, the high injection of electron hole pairs may also have an effect. For this condition, the peak concent ration needs to be more than 1019 cm 3 (as shown by Figure 5 6) in order for the piezoresistance factor to be reduced. Even then, the immediate drift and eventually diffusion following the particle strike will quickly reduce the peak concentration carriers and likely minimize any effects on the piezoresistance factor. This issue requires further experimental investigation for which carrier concentrations above 10 19 cm 3 5.3.4 Piezoresistive Mobility Model Implementation As discussed in the previous section, pi ezoresistance defines the relationship between electric field, current, and mechanical stress The relationship between piezoresistance and mobility formulated as PAGE 131 131 ( 1 86 ) or in terms of current density J as ( 1 87 ) where is resistivity, is piezoresistance, is stress (not to be confused with conductivity) and is mobility. E quation ( 1 86 ) assumes that mobility changes li nearly with stress. This relationship is reasonable as long as stress value remains below ~1 GPA since experimental data show a linear trend for mobility versus stress for both n type and p type silicon [ 67 ]. Since stress does not typically exceed ~1 GPA f or a MOSFET channel (45 nm node), this is a reasonable assumption to make for this work. However, non linear piezoresistive modeling will need to be considered for future technology nodes. For the general case of stress in an unprimed coordinate system, the change in mobility due to stress (less than ~1 GPa) is ( 1 88 ) or in tr ansformable (orientation) coordinate system ( 1 89 ) PAGE 132 132 The transformed piezo orientation of [110] on a (001) wafer for bulk silicon are given below. For electron mobility, it can be written as ( 1 90 ) For hole mobility, the matrix can be written as ( 1 91 ) For clarity, the subscripts in the two previous matrixes are equivalent to (1=X=[110]), (2=Y=[1 10 ]), (3=Z=[001]), and so on. As discussed earlier, the full set of piezoresistance coefficients as a function of orientation have been derived in the Appendix such that equation ( 1 89 ) can be used for any silicon orientation. The change in current density due to stress can be written as ( 1 92 ) where t he above equation is reduced to the following PAGE 133 133 ( 1 93 ) Expanding each current component, equation ( 1 93 ) reduces to ( 1 94 ) ( 1 95 ) ( 1 96 ) 5.4 Uniaxial Strained S i Diode Although strained Si technology has been widely adopted, the effects of mechanical stress on current transients generated by laser or ion strikes at the source/drain regions had not been reported until recently by Park, Cummings, Arora, and colleag ues [ 42 ]. It is important to understand how mechanical stress affects these transient pulses since the transport of the radiation generated carriers in the substrate is affected by stress. Laser induced current transients on a uniaxially stressed Si N+/P junction diode are discussed in this section [ 42 ]. An N+/P diode is a good representation of the source/drain junctions that are responsible for charge collection in n channel MOSFETs. Furthermore, stress induced electron mobility enhancement is easier to understand than that of holes [ 47 ], so N+/P diodes were used in this work. P channel MOSFETs are also important for considering single event transients but will be discussed more in the next section. The shapes of current transients and the amount of coll ected charges are measured as a function of stress, because both of them are PAGE 134 134 crucial in predicting SEUs in circuits [ 1 ]. Additionally, the results of the diode experiments and simulations will give some insight in how strained silicon technology affects si ngle event transients in modern CMOS devices. 5.4.1 Experimental Setup Controlled external mechanical stress was applied via a four point bending jig [ 56 ] while the samples were irradiated using a picosecond pulsed laser as in Figure 5 7. The samples used in th is study are N+/P diodes fabricated on (001) Si wafers using a standard 130 nm CMOS technology. The active area of the diodes is 50 m 100 m. Nickel silicide (NiSi), silicon oxide (SiO x ), and copper (Cu) patterns are present on top of the diodes as show n in Figure 5 7 using transmission electron microscopy (TEM) and energy dispersive X ray spectroscopy (EDS). The thickness of the NiSi, SiO x and Cu patterns is ~20 nm, 720 nm, and 280 nm, respectively. The doping densities of the n+, p well, and p substra te are ~10 20 ~10 18 and ~10 16 cm 3 respectively [ 57 ]. A cavity dumped dye laser with a wavelength of 590 nm, a pulse energy of 218 pJ, and a pulse width of 1 ps is used to inject electron hole pairs in the diode. The laser direction is normally incident to the diode surface and has a spot size of 12 m in diameter. The peak carrier concentration produced by the laser is ~1.6 10 19 cm 3 The pulse laser energy reaching the diode active area is smaller than the value measured at the surface of the structur e due to the optical properties of the layers on top of the diode [ 58 ] [ 59 ]. Current transients on the N+/P diode are measured under different values of stress (160 MPa and 240 MPa tensile, no stress, and 160 MPa compressive) with a 5 V reverse bias. The experimental setup and analysis are discussed in greater deal in [ 42 ]. PAGE 135 135 5.4.2 Comparison of Experimental and Simulation Results The FLOODS simulation tool was used to explain the mechanisms responsible for the differences in charge collection between stresse d and unstressed devices. Additionally, the simulations were used to predict the effects of high mechanical stress (~1 GPa) on laser induced current transients, above the maximum stress that could be applied using the four point bending jig (240 MPa). Ba sed on the experimental analysis discussed in the previous section, FLOODS simulations were performed to understand the mechanisms of carrier transport under uniaxial stress and to predict how high stress (~1GPa) affects the current transients in diodes. T he Masetti and Brooks Herring mobility models were used to account for carrier transport in a high injection case (note: the high injection mobility model in Chapter 6 was not available at the time of this study). Shockley Read Hall and Auger band to band recombination models were also considered. The number and distribution of electron hole pairs generated by the laser pulse was calculated by a single photon absorption (SPA) equation discussed in chapter 2 [ 22 ]. The c hange in the amount of generated e h ca rriers was calculated to be less than 3% for 1 GPa of uniaxial tensile stress [ 42 ]. Before analyzing the effects of stress on current transients, baseline simulations under no stress were performed. These results were matched to the measured current transi ent under no stress. It is very important to understand the physics that dominates current transients in an unstressed case in order to predict the results under a stressed case. A 2 dimensional simulation structure, shown in Figure 5 9, was built based on analysis of the structure and material of the N+/P diodes, as discussed in [ 42 ]. The width and depth of the diodes were set to 100 m and 10 m, respectively to prevent carrier reflection at the boundaries. The piezoresistive mobility model (discussed ear lier PAGE 136 136 in this chapter) coefficients was used to consider mobility enhancement under mechanical stress [ 60 ]. Additionally, t he doping dependence of the coefficients is considered [ 61 ] The simulated current transien ts in Figure 5 10 show the same trend as the experimental d ata in Figure 5 11 I max and Q in the simulations also agreed with the experiments, as shown in Figure 5 12 and Figure 5 13 The data points in the experiments are the average I max and Q at each stress level. The error bars in the data points represent the standard deviation in the data at each stress level. The simulation results predicted that I max and Q under 1 GPa of tensile stress will decrease by ~23% and ~21%, respectively. Analogous to tensile stress, 1 GPa of compressive stress increased I max and Q by 17% and 13%, respectively. The experiment and simulation results for strained N+/P diodes showed that uniaxial stress changes the shape of current transients and collected charges. 5.4.3 Uniaxially Strained Si Diode Summary This section showed that uniaxial tensile stress in Si N+/P diodes decrease the maximum peak currents and collected charges for laser induced current transients. Quantitative analysis and FLOODS simulation results suggest that this can be attributed to the degradation of electron mobility along the [ 001 ] direction. In other words, the zz zz component in equation ( 1 89 ) Using Figure 5 14 as a reference, consider the stress and piezoresistance zz zz zz zz component can be expanded into the following term as ( 1 97 ) PAGE 137 137 Due to the orientation of the diode, the zz term is equivalent to 33 Thus, the previous equation can be written in terms of the x y and z axis as ( 1 98 ) Because the yy and zz stress components are negligible when uniaxia l stress is point bending jig, equation ( 1 98 ) can be reduced to the following as ( 1 99 ) where the subscript n represents electron mobility. This shows that a positive stress (tensile) for xx will reduce the electron mobility whereas a negative stress (compressive) will increase the mobility in the [001] direction. Additionally, this result is verified by Figu re 5 11 Therefore, uniaxial strain engineering has the potential to control the shape of single event transients and the amount of charges collected in devices. This will be explored more in the next section for CMOS devices. 5.5 Predictions for Strained Si MOSFETs The results of the uniaxially strained diode in the previous section can be extended to the modern CMOS technology. Uniaxial strained silicon is considered in this section since it is a leading technology for enhancing transistor performance for su b 100 nm logic technology [ 62 ], [ 63 ]. Additionally, uniaxial mechanical stress improves device characteristics such as mobility and gate tunneling current, with minimal stress induced threshold voltage shifts [ 64 ]. Building upon the N+/P diode work, this s ection investigate s how strained Si technology impacts charge collection and current transients for 45 nm CMOS devices. PAGE 138 138 5.5.1 Simulation Setup Overview It is necessary to perform both process and devices simulation for this section. Front end process simulations are required to calculate the stress contours in all directions. Unlike the uniaxially strained diode (stress only along <110> direction), high stress values often occur in every direction for a modern CMOS device. The main focus of the FLOOPS process sim ulations was to closely model TSMC production level CMOS process at the 45 nm node [ 48 ]. These devices were modeled since data on the process, structural dimensions, and current voltage characteristics were readily available [ 48 ]. For 45 nm node devices, the CMOS channels are oriented in the <110> direction since it is advantageous for mobility enhancement. In order to induce advantageous stress along the channel, a tensile capping layer is used for the NMOS devices and embedded SiGe with a compressive cap ping layer is used for PMOS devices, as shown by the schematic in Figure 5 15 [ 48 ]. Because g ermanium is larger than silicon when it sits on a substitution al lattice site a local lattice expansion occurs [ 53 ]. At high concentrations significant strain va lues can result due to a lattice mismatch between the silicon substrate and the dopants as shown in Figure 5 16 This approach to inducing stress in the channel is very common and is discussed in other work [ 65 ], [ 66 ]. In addition to strained silicon proce sses, the shallow trench isolation (STI) regions are designed to have low stress and the gate equivalent oxide thickness (EOT) is about 1.5 nm. A TEM example of the CMOS 45 nm node from other is given in Figure 5 17 where it can been seen that the typical gate length is about 30 nm [ 66 ],[ 48 ]. A fully processed 2 D MOSFET is shown in Figure 5 18 The gate, oxide, spacers and capping layer processes (deposition, etching, etc.) were simulated by FLOODS. These geometries are all factors in how stress is calcul ated for the NMOS and PMOS PAGE 139 139 devices. Typically, stresses of around 1 GPa have been reported in the channel for the 45 nm node and the stress inducing processes for the simulations were designed to induce such stress. Figure 5 18 also shows the boundary of t he device which is 0.8 5 length (for 3 D simulations). More importantly, the structure is large enough to bound the entire strike path. For the particle strike, a 1 MeV h elium ion (a.k.a. alpha particle) is used to generate the single event transient and uses the Gaussian profile given by equation ( 1 3 ) At this energy, the io n has a stopping range of 3.54 micrometers and an LET of 1.312 MeV cm 2 /mg, as calculated by SRIM. The alpha particle is a useful illustration since these particles are becoming increasingly problematic as devices are downscaled [ 1 ]. To create a worst case scenario, the particle strike is path was setup to go directly through the drain region, about 150 nm away from the center of the gate. For a clearer visualization, the 3 D strike path is shown in Figure 5 19 where the 10 18 cm 3 charge contour is shown. F or the device simulations, the general purpose mobility model discussed in Chapter 6 was used. Additionally, velocity saturation (Canali model) and transverse gate field effects (Lombardi model) were included. For recombination, the Auger and SRH models we re used. The quasi Fermi discretization approach was used since the stress calculations in FLOOPS are performed using a finite element approach. Diffusive boundary sinks were used on the device edges to minimize carrier reflection. The devices were biased to V DS =1.0 V (NMOS) and V DS = 1.0 V (PMOS). PAGE 140 140 Prior to simulating the short channel MOSFETS (L g =30 nm), long channel devices (L g enhancement). A uniaxial stress of 1 GPa was induced in the <110> direction and the linear current enhancement ( D,LIN /I D,LIN,0 ) was found to be ~32% for the NMOS and ~73% for the PMOS. This result agrees well with the piezoresistance coefficients for both the NMOS ( 31.210 11 /Pa) and PMOS (71.810 11 /Pa) devices in the <110> direction. The measured current voltage characteristics for the TSMC 45 nm CMOS devices are shown in Figure 5 20 [ 48 ]. FLOODS device simulations (including process induced stress) were performed for the previously described CMOS devices in Figure 5 18 The current voltage characteristics for the NMOS devices are shown in Figure 5 21 and Figure 5 22 where it can be seen that the results agree very closely with the experimentally measure devices. For a tensile channel stress of ~1 GPa, the NMOS drain current enhancement is about 14%. Next, the current voltage characteristics for the PMOS devices are shown i n Figure 5 23 and Figure 5 24 For a compressive channel stress of ~1 GPa, the PMOS drain current enhancement is about 19%. It should be noted that the enhancement for short channel devices is lower than long channel devices as shown by Figure 5 25 [67 ]. The physical mechanisms for this behavior are still under investigation. 5.5.2 NMOS Simulation Results Before discussing the simulation results, it will be useful to have a visual reference for the MOSFET orientation. Figure 5 26 shows the MOSFET orientation whe re the channel is aligned in the <110> direction. For example, the stress xx component is in PAGE 141 141 zz zz in the direction of the charge strike, the zz component is in the The FLOOP S predicted stress profiles for the 2 D NMOS simulations are shown in Figure 5 27 and Figure 5 28. Although the tensile capping layer induces a significant in the For the strike region, both the xx and zz are quite small zz zz contribution quite small for the NMOS. The contributions to the charge strike (in 2 zz zz direction <001> are given by ( 1 100 ) which can be derived from equation ( 1 89 ) As with the uniaxially strained diode, the change in mobility in [001] direction has the largest impact on charge collection and current; the single event current flow is primary in this direction due to the depletion region and funneling field. Very little change in the current transient and charge collection are observed for the strained silicon NMOS as shown by Figure 5 29 and Figure 5 30 x x zz yy [1 10 ] component should also be considered since it contributes to the zz mobility in the direction of the strike. The FLOOPS predicted stress profiles for the 3 D NMOS simulations are shown in Figure 5 31, Figure 5 32 and Figure 5 33. Altho ugh the tensile capping layer induces a significant amou nt of tensile stress in the channel direction, only a fraction of the xx yy and zz components zz zz component PAGE 142 142 small as we ll for the NMOS. The contributions to the charge strike (in 3 D) for the zz zz direction are given by ( 1 101 ) which can be derived from equation ( 1 89 ) Very little change in the current transient and charge collection are observ ed for the strained silicon NMOS as shown by Figure 5 34 and Figure 5 35 5.5.3 PMOS Simulation Results The FLOOPS predicted stress profiles for the 2 D PMOS simulations are shown in Figure 5 36 and Figure 5 37. The compressive capping layer and embedded SiGe in duces a significant amount of compressive stress in both the channel direction and the depth direction. For the strike region, both the xx and zz are large near the drain zz zz component significance for the upper portion of th e strike path. The contributions to the charge strike (in 2 zz zz direction are given by ( 1 102 ) which can be derived from equation ( 1 89 ) A slight increase in the current transient peak and charge collection are obs erved for the strained silicon PMOS as shown by Figure 5 38 and Figure 5 39. However, this increase does not include the yy component of stress since the results are for a 2 D simulation. In addition to the x [110] and z [001] directions, the y [1 10 ] sho uld also be zz mobility component. The FLOOPS predicted stress profiles for the 3 D PMOS simulations are shown in Figure 5 40, Figure 5 41 and PAGE 143 143 Figure 5 42. The compressive capping layer and embedded SiGe induces a s ignificant amount of compressive stress in both the channel, perpendicular and the depth directions. For the strike region, both the xx yy and zz are significant near the drain junction. The contributions to the charge strike (in 3 zz zz direction are given by ( 1 103 ) which can be derived from equation ( 1 89 ) Very little change in the current transient and charge collection are observed for the strained silicon PMOS as shown by Figure 5 43 and Figure 5 44. This is due to the fact that for the [110] channel orientation, the zz zz = = 1.1 =6.6 [10 11 Pa]). Because the piezoresistance coefficients are so small, high stress values will do zz zz component. Thus, stress will always have a minimal impact on charge collection for PMOS device on a (001) wafer oriented in the [110] direction. 5.5.4 Impact of STI on Single event Transients In the previous section, it was shown that strained silicon has only a minor impact on single event behavior for typical CMOS devices at the 45 nm node. This was mainly due to the fact that the process induced stress was isolated near the surface in the channel and source/drain regions. In con trast, the uniaxially strained N+/P results show that stress has a large impact on single event transients (SET) since the stress profile goes deeper into the bulk. A deeper stress profile can result in a more significant change in SET results. PAGE 144 144 One possib le way to induce stress deeper into the substrate for modern CMOS devices is to use shallow trench isolation (STI) techniques. During front end processing, stress is generated between the STI regions (i.e. source, drain, channel) due to the lattice mismatc h created at the STI silicon sidewall interface. Much research has been performed to understand the effect that STI has on mobility, saturation velocity, and threshold voltage [ 67 ] For typical 45 nm CMOS fabrication, it was found that STI induced disadvan tageous stress (less ideal mobility enhancement) in the channel region [ 68 ]. For example, any compressive stress (resulting from STI) along the channel of a NMOS device will reduce the mobility. Thus, modified STI processing techniques were developed to mi nimize the stress in the STI regions [ 69 ]. Alternative to these techniques, some approaches considered using STI regions to induce advantageous stress in the channel in lieu of embedded SiGe and capping layers [ 70 ], [ 71 ], [ 72 ] This has strong implications for SET behavior because STI processes can induce stress deeper into the substrate than capping layers. In this section, 3 D simulations are run to compare the effect of strained and unstrained (or minimized) STI regions on charge collection. 5.5.5 STI Simulat ion Results Similar to the previous process simulations for the CMOS devices with SiGe and capping layers, a lattice mismatch can be created in the STI regions to induce stress. Work by several research groups has shown that both compressive and tensile st ress could be induced with STI using a high aspect ratio process (HARP) with a O 3 / tetraethoxylonesilane (TEOS) based subatmospheric chemical vapor deposition (SACVD) trench fill process [ 70 ],[ 72 ]. Figure 5 45 shows the possible stress values that can be in duced using the aforementioned processing techniques [ 70 ]. Interestingly, it is PAGE 145 145 shown that tensile stress induced from STI regions can theoretically reach up to 1 GPa for the NMOS channel PMOS results are not shown since the piezoresistance coefficients a re very small for the <001> direction. A three dimensional NMOS simulation structure, identical to the one in the previous section, was used for the processing, steady state device and transient simulations. The depth of STI was chosen to be 350 nm in or der to match the 2007 ITRS guidelines for the 45 nm node [ 39 ]. The structural layout, size, and particle strike model (alpha particle) are the same as the previous MOSFET simulation section. The 3 D stress profiles generated by FLOOPS are shown in Figure 5 46, Figure 5 47, and Figure 5 48. For the sake of argument and comparison, the simulations assume that up to 1 GPa of stress can be induced between STI regions. As expected, the stress profiles go much deeper into the device and these results agree with o ther work [ 70 ]. Interestingly, the stress components in every direction work to reduce the electron mobility in the [001] direction. This results in a significant reduction of the current peak and collected charge as shown by Figure 5 49 and Figure 5 50 To gain insight into the simulation results for all the above mentioned strained silicon devices, consider Figure 5 51 For the case of electron mobility, the tensile capping layer only induced stress at the surface so the resulting change in mobility in [ 001] was minimal. However, a larger result was seen for the uniaxially strained diode because the stress profile was uniformly deep into the substrate. Finally, the largest change in electron mobility is seen for the STI induced stressed NMOS. A high amoun t of stress went deeper into the device and the resulting change in electron mobility was large. For the STI case, the depth of the stress profile was deep enough to bound the PAGE 146 146 funneling region and thus the single event results are significantly different t han those for an unstrained device. For the PMOS devices, the piezoresistance coefficients in the <001> are too small to induce any sort of meaningful change in hole mobility as shown in Figure 5 52. However other orientations should be explored for PMOS devices in order to exploit higher piezoresistance values that are in other directions. 5.5.6 Strained Si MOSFET Summary The results for the strained silicon NMOS and PMOS devices show that the depth of the stress profile is very important for single event effe cts. For the devices that used SiGe and capping layers to induce stress, the change in charge collection was minimal since the stress was limited to the surface. However, for an NMOS with STI induced stress, the stress profile was much deeper into the subs trate. Since particle strike paths can go deep into the bulk of a device, a deeper stress profile (thus mobility change) will have a larger impact on collected charge. Predictive simulation results for an NMOS with 1 GPa of STI induced stress show that a ~ 30% reduction in charge collection and current can be attained. Such knowledge can be useful for mitigating the effects of SEU for modern devices. T he results sugg est that strained Si technology could have a significant impact on SEUs at the circuit level. PAGE 147 147 Table 5 1. Values of piezoresistance ( ) coefficients ( 10 11 Pa 1 ) used in FLOOD S [60] Si 0 11 12 44 n type 11.7 102.2 53.4 13.6 p type 7.8 6.6 1.1 138.1 PAGE 148 148 Figure 5 1 Ellipsoids of co corresponding to one of the degenerate conduction band valleys. A) Unstrained Si. B) S trained Si C) Energy level at the bottom of the six conduction band valleys. Advantageous strain splits the energy levels as shown [ 51 ] Figure 5 2 Simp lified schematic of strain induced hole energy band splitting and the intervalley phonon scattering process. [ 51 ] PAGE 149 149 Figure 5 3 Linear spring element in equilibrium (top) and then subjected to tensile forces (bottom). Figure 5 4 Three dimensional stresses on an element. Figure 5 5 Baseline t e nsor orientation notation (and Miller indices) for silicon. [110] 45 PAGE 150 15 0 Figure 5 6 Piezoresistance factor P ( N,T ) as a function of impurity concentration and temperature f or n type Si [ 61 ]. Figure 5 7 Laser induced current transient measurement system using a four point bending jig. [ 42 ] PAGE 151 151 Figure 5 8 Schematic of N +/P diode structure throu gh TEM and EDS analysis ( top ) and TEM image (bottom) [ 42 ] PAGE 152 152 Figure 5 9 Schematic of laser induced current transients and 2 dimensional simulation structure of an n+p diode. [ 42 ] Figure 5 10 Laser induced current transients and the ratio of collected charge measured as a function of < 110> uniaxial mechanical stress [ 42 ]. PAGE 153 153 Figure 5 11 Simulated laser induced current transients as a function of <110> uniaxial mechanical stress [ 42 ] Figure 5 12 Peak current (Imax) as a function of mechanical stress. (positive (+) : tensile, negative ( ): compressive) [ 42 ]. PAGE 154 154 Figure 5 13 Collected charges (Q) until 10 ns. (positive (+) : tensile, negative ( ): compressive) [ 42 ]. Figure 5 14 Orientation for the N+/P diode exper iment and simulations. PAGE 155 155 Figure 5 15 Strained Si CMOS technology for 45 nm node. CESL represents the compressive (PMOS) and tensile (NMOS) capping layers [ 48 ]. Figure 5 16 Lattice expansion from germanium [ 53 ]. Figure 5 17 TEM micrographs of 45 nm node transistors. A) NMOS [ 65 ]. B) PMOS [ 65 ]. C) e SiGe PMOS [ 48 ]. PAGE 156 156 A B Figure 5 18 2 D simulation structure. A) 2 D MOS device after processing in FLOOPS. B) MOS device boundary and strike path. Boundary sinks (discussed in Chapter 4) were used on the right and left ( and front and back for 3 D) device edges. A B Figure 5 19 3 D MOSFET structure and Helium particle strike path. The 10 18 cm 3 charge contour is shown in green. A) 3 D mesh. B) 3 D particle strike PAGE 157 157 A B Figure 5 20 Measured I V characteristics for 45 nm strained Si CMOS. A) I D V GS characteristic. B) I D V DS characteristic [ 48 ]. Figure 5 21 FLOODS predicted I D V GS characteristic for a strained si licon NMOS device (45 nm). PAGE 158 158 Figure 5 22 FLOODS predicted I D V D S characteristic comparing a strained and unstrained NMOS device (45 nm). I D SAT enhancement is about 14% (~1 GPa tensile channel stress). Figure 5 23 FLOODS predicted I D V GS characteristic for a strained silicon PMOS device (45 nm). PAGE 159 159 Figure 5 24 FLOODS predicted I D V D S characteristic comparing a strained and unstrained PMOS device (45 nm). I D SAT enhancement is about 19% (~1 GPa compressive channel stress). Figure 5 25 I D,lin enhancement versus uniaxial longitudinal tensile stress plotted fo r 10 and 0.1 67 ] PAGE 160 160 Figure 5 26 MOSFET orientation (and associated notation) with the channel in the <110> direction. Figure 5 27 NMOS Str ess XX component (channel direction) in [Pa] units. 2 D FLOOPS simulation results. A tensile capping layer induces a tensile stress (~ 1 GPa) in the NMOS channel. Strike path shown by arrow. PAGE 161 161 Figure 5 28 NMOS Stress ZZ component (depth direction) in [ Pa ] units 2 D FLOOPS simulation results. A tensile capping layer induces very little stress in the depth direction <001> Strike path shown by arrow. Figure 5 29 2 D NMOS current transient for strained and unstrained devices. V ds =1.0 V. PAGE 162 162 Figure 5 30 2 D NMOS charge collection for strained and unstrained devices. V ds =1.0 V. Figure 5 31 NMOS Stress XX component (channel direction) in [Pa] units. 3 D FLOOPS simulation results. A tensile capping layer induces a tensile stress (up to 1 GPa) in the NMOS channel. Strike path shown by arrow. PAGE 163 163 Figure 5 32 NMOS Stress YY component ( perpendicular to channel) in [Pa] units. 3 D FLOOPS simulation results. A tensile capping layer induces lower stress (~100 500 MPa) perpendicular to the NMOS channel. Figure 5 33 NMOS Stress ZZ component (depth direction) in [Pa] units. 3 D FLOOPS simulation results. A tensile capping layer induces very little stress in the depth direction <001>. PAGE 164 164 Figure 5 34 3 D NMOS current transient for strained and unstrained devices. V ds =1.0 V. Figure 5 35 3 D NMOS charge collection for strained and unstrained devices. V ds =1.0 V. PAGE 165 165 Figure 5 36 PMOS Stress XX component (channel direction) in [Pa] units. 2 D FLOOPS simulation results. A compressive capping and embedded SiGe layer induces a compressive stress (up to 1 GPa) in the PMOS channel. Strike pa th shown by arrow. Figure 5 37 P MOS Stress ZZ component (depth direction) in [ Pa ] units 2 D FLOOPS simulation results. A compressive capping layer and embedded SiGe induces significant compressive stre ss in the depth direction <001> PAGE 166 166 Figure 5 38 2 D PMOS current transient for strained and unstrained devices. V ds = 1.0 V. Figure 5 39 2 D PMOS charge col lection for strained and unstrained devices. V ds = 1.0 V. PAGE 167 167 Figure 5 40 P MOS Stress XX component (channel direction) in [Pa] units. 3 D FLOOPS simulation results. A compressive capping layer and embedded Si Ge induces a compressive stress (up to 1 GPa) in the P MOS channel. Strike path shown by arrow. Figure 5 41 P MOS Stress YY component ( perpendicular to channel) in [Pa] units. 3 D FLOOPS simulation results The compressive capping layer and embedded SiGe induces lower stress (~1 GPa) perpendicular to the P MOS channel. Strike path shown by arrow. PAGE 168 168 Figure 5 42 P MOS Stress ZZ component (depth direction) in [P a] units. 3 D FLOOPS simulation results. A compressive capping layer and embedded SiGe induces significant compressive stress in the depth direction <001>. Strike path shown by arrow. Figure 5 43 3 D P MOS current transient for strained and unstrained devices. V DS = 1.0 V. PAGE 169 169 Figure 5 44 3 D P MOS charge collection for strained and unstrained devices. V ds = 1.0 V. Figure 5 45 Hysteresis effect of the deposited film as a function of temperature (nitrogen ambient). The stress of the film is fully stable after the first anneal cycle. [ 70 ] PAGE 170 170 Figure 5 46 N MOS Stress X X component (channel direction) for STI induced stress. Figure 5 47 N MOS Stress YY component ( perpendicular direction) for STI induced stress. 3 D FLOOPS simulation results. PAGE 171 171 Figure 5 48 N MOS Stress ZZ component ( depth direction) for STI induced stress. 3 D FLOOPS simulation results Figure 5 49 3 D N MOS current transient for STI strained and unstrained devices. V DS =1.0 V. PAGE 172 172 Figure 5 50 3 D N MOS collected charge for STI strained and unstrained devices. V DS =1.0 V. Figure 5 51 Electron mobility change along the particle strike path in the <001> direction as a function of depth for the 3 D NMOS device. PAGE 173 173 Figure 5 52 Hole mobility change along the particle strike path in the <001> direction as a function of de pth for the 3 D PMOS device. PAGE 174 174 CHAPTER 6 BULK MOBILITY MODELI NG FOR SINGLE EVENT EFFECTS 6.1 Introduction Mobility is a key parameter in characterizing electron and hole transport in semiconductor devices. The results of semiconductor device simulations are highly depe ndent on the accuracy of the mobility models used. For instance, the overall effect of mobility on current density can be shown in terms of quasi Fermi levels as ( 1 104 ) ( 1 105 ) where n and p are the electron and hole densities, n,p the quasi Fermi levels, J n,p the current density and n,p the mobilities. Therefore, it is important to choose an accurate mobility model so that the simulation results will be relevant. Mobility in silicon is controlled by scattering, it is commonly expressed as ( 1 106 ) where m is the mean free time between collisions and m* the conductivity effective mass [ 49 ]. Because there are multiple scattering mechanisms in silicon (i.e., ionized m can be defined in terms of the individual mean free times by ( 1 107 ) Since mobility is proportional to the mean free time as in equation ( 1 106 ) it can be formulated in terms of each of these scattering mechanisms. By using the Mattheissen rule and following the same form as ( 1 107 ) bulk silicon mobility can be formulated as PAGE 175 175 ( 1 108 ) where the different components of the mobility are represented by i and the total effective mobility is T The most significant bulk silicon mobility contributions are from scattering from the lattice, donor, acceptor and carrier carrier inte ractions Although there are many different approaches to modeling mobility in silicon, most models use the form of equation ( 1 108 ) to account for all the sc attering mechanisms since the equation is computationally efficient and reasonably accurate. However, most models only account for a few mechanisms at a time. Therefore, it is desirable to combine the most accurate dependencies (e.g., doping levels, temper ature, carrier carrier scattering) from existing mobility models into a single mobility model set suitable for device simulations. The manner in which mobility at high injection levels is modeled is especially important since a large number of electron hol e pairs are generated along a particle strike path. Since a particle strike generates an equal number of free holes and electrons, the mobility is qualitatively important because it affects how rapidly and how far the deposited charges separate, and hence has a first order effect on the potential distribution and charge collection during the strike recovery. Chapter 2 gave a brief example of the impact of mobility on the total charge collection and transient current characteristic. In this chapter, mobility will be discussed in much greater detail where the focus will be on modeling mobility in the bulk region of the device since that area is important for SEE. This chapter starts by giving an overview of existing mobility models commonly used for device si mulations. Next, two proposed mobility models are formulated and PAGE 176 176 tested. Each model is in a computationally efficient form and accounts for majority and minority carrier mobility, carrier carrier scattering and temperature dependence. Finally, several fiel d dependent models important for CMOS simulations are discussed. These models account for lateral (channel direction) field effects such as velocity saturation and transverse field effects such as surface roughness and surface phonon scattering. 6.2 Overview of Existing Bulk Mobility Models Due to the large number of free carriers that exist in the substrate immediately following a particle strike, it is important to model carrier mobility in the bulk of the device. For radiation effects simulations, various bulk mobility models for device simulation are available. A thorough summary of conventional mobility models is given in Figure 6 1 which shows that a wide variety of models are available for bulk silicon, each with particular advantages for device simula tion. Some models focus on the accurate fitting of majority mobility versus doping levels, some on minority mobility and others on temperature dependence. Each model is qualitatively compared against others with respect to majority carrier mobility, minori ty carrier mobility, electron hole scattering, screening of charge carriers, and temperature dependence. As evident in Figure 6 1, no single model accurately accounts for every parameter. For example, the Masetti model can be used for its excellent fitting to majority carrier data but lacks a carrier carrier scattering description, limiting its applicability in situations with high carrier densities, e.g., following an ion strike [ 72 ]. Furthermore, very few models focus on the electron hole scattering mecha nism, which is important for simulating radiation effects, such as single event upsets. An important aspect of radiation effects simulations is how the mobility model treats high injection electron hole carrier densities. As pointed out by Dodd [ 15 ], the PAGE 177 177 charge densities immediately after the passage of an ionizing particle can exceed 10 20 cm 3 For carrier densities below 10 18 cm 3 Dannhauser [ 73 ] and Krausse [ 74 ] measured the sum of electron and hole mobilities as a function of the concentration of carr iers injected into the weakly doped region of a silicon P I N diode. Unfortunately, very little experimental data has been published for electron hole carrier densities above 10 1 8 cm 3 Although limited data are available, approximations based on semi clas sical quantum theory, such as the Conwell Weisskopf or Brooks Herring models, predict that an increase in electron and hole density results in a decrease in carrier mobility [ 75 ]. Two bulk mobility models that account for carrier carrier scattering are the Philips unified mobility model and the Dorkel Leturcq mobility model. The Philips unified mobility model is a commonly used mobility model for device simulation and has been used for recent simulation work in the area of CMOS and SiGe HBT radiation effec ts [ 76 ],[ 77 ]. T he Philips model accounts for majority and minority carrier mobility, the screening of the impurities by charge carriers, electron hole scattering, clustering of impurities, and temperature dependence [ 37 ] However, t he carrier carrier scatt ering in the Philips model is formulated in such a way such that it does not match known experimental data for electron and hole concentrations above 10 1 7 cm 3 Therefore, TCAD simulations result in single event current pulses that are too large when using the Philips model, and hence voltage pulse widths that are too short as discussed in [ 15 ] For single event simulations, t he Dorkel Leturcq model has been suggested as a better alternative to the Philips model since at high electron hole densities, the mo bility agrees better with measured data [ 4 ] This model describes mobility in terms of doping PAGE 178 178 dependence and carrier carrier scattering. However, for modern devices it lacks accurate majority and minority mobility descriptions since the model was primarily designed for doping levels below 10 19 cm 3 [ 78 ] A lso, a disadvantage of the Dorkel Leturcq model is that it does not fit the data well at high doping concentrations and has not been formulated for minority carrier mobility Due to inconsistencies betwee n existing bulk models and experimental data, alternative approaches to modeling mobility are presented in the next two sections. The proposed models account for majority and minority carrier mobility and temperature dependence in a computationally efficie nt form. First, a high injection mobility model (a.k.a. UF mobility model) is formulated to specifically to account for electron hole scattering that occurs during a particle strike. Next, a general purpose model is formulated to address some of shortcomin gs of the UF mobility model and to account for the screening of charge carriers. 6.3 High Injection Mobility Model The goal of the high injection mobility model is to formulate a mobility model suitable for radiation effects simulations that accurately descri bes majority and minority carrier mobilities, carrier carrier scattering, and temperature dependences [Cum10a] There are several ways to approach the modeling of mobility. Some methods formulate mobility starting from fundamental quantum mechanics princip les and therefore are very computationally intensive [ 80 ] Other mobility modeling methods start with simplified formulations of lattice and ionized impurity scattering (as discussed previously) and then use fitting paramet ers to match experimental data T he UF mobility model use s the latter approach to modeling mobility since computational efficiency is important for device simulations. As discussed in the following sections, the proposed model PAGE 179 179 combines the most accurate dependencies (e.g. doping levels, temperature, carrier carrier scattering) from existing mobility models to form a single mobility model set suitable for radiation effects device simulations in silicon. 6.3.1 Majority Carrier Modeling The majority carrier modeling in this section describes the lattice scattering and ionized impurity scattering processes of electrons in n type material and holes in p type material. To formulate the majority carrier mobility for the proposed model, the well defined doping and temperature functions in the Masetti a nd Arora models will be combined. The mobility derivation is best understood by starting with the modeling approach of Caughey Thomas which shows that plots of experimentally measured mobility data versus the logarithm of doping density strongly resemble t he Fermi Dirac function [ 84 ]. The Caughey Thomas mobility model in terms of doping density is expressed as ( 1 109 ) where C ref and are fitting parameters, N is the total doping density, and min and max model is suitable for lower impurity concentrations but is inaccurate at higher concentrations. Building upon ( 1 109 ) a third term is added to account for the additional decrease in mobility that occurs when the doping level is more than 510 19 cm 3 [ 72 ]. This results in the Mas etti mobility model and is of the form ( 1 110 ) PAGE 180 180 The Masetti model is s hown in Figure 6 2 and Figure 6 3 where it is compared against the Dorkel Leturcq, Philips and proposed mobility models. The Masetti model has been fitted to experimental data very accurately for both electrons and holes since majority carrier mobility has been heavily investigated. The parameters for the majority carrier mobility are given in Tabl e 6 1 and are based on [ 72 ]. A disadvantage of the Masetti formulation is that it is not a function of temperature. To add temperature dependence, the Arora mobi lity model approach is used since it is well fit to experimental data with mobility as a function of temperature [ 81 ]. The Arora model can be formulated in terms similar to the Caughey Thomas expression in ( 1 109 ) where the terms min max C ref and can be written as functions of temperature [ 81 ],[ 21 ]. Using the same approach, but building on the Masetti formulation in ( 1 110 ) the new proposed majority carrier mobility can be written as ( 1 111 ) where T n =( T /300 K). The sub script i stands for e (electrons) or h (holes) and the T n terms are the temperature fitting parameters. The third term on the right hand side of ( 1 111 ) is not a function of temperature since for high impurity conc entrations, the carrier mobility in silicon becomes nearly temperature independent [ 82 ]. The values for the temperature fitting parameters are given in Table 6 2. The parameters are based on Arora's model but are modified to fit the experimental temperatur e data in [ 82 ],[ 83 ],[ 85 ] A comparison between the proposed model, the Arora model, and measured data for both electron and hole mobilities is given in Figure 6 4 and Figure 6 5. The plots PAGE 181 181 show that both models follow a similar mobility trend over a range of temperatures. Since the Arora model uses a formulation based on Li and Thurber [ 82 ] and the proposed model follows Masetti [ 72 ], a small difference in results is observed. For doping levels higher than 10 19 cm 3 the proposed model fits experimental dat a better since the Arora model over predicts mobility at high doping levels. 6.3.2 Minority Carrier Modeling Minority carrier mobility is a description of the scattering processes of electrons in p type material and holes in n type material. As with the majority carrier formulation in the previous section, a similar approach is used to model the minority carrier mobility by using the Caughey Thomas and Masetti expressions as a starting point. Because the Masetti model does not include minority carrier mobility, a new set of fitting parameters is used. Following th e temperature dependence approach in ( 1 111 ) the new proposed formulation for minority carrier mobility i s of the form ( 1 112 ) where fit is an additional fitting term. The fo urth term on the right hand side of ( 1 112 ) arises from the fact that experimental data show that minority carrier mobility exceeds majority carrier mobility at high doping concentrations (~110 18 110 20 cm 3 ) [ 86 ]. This additional fitting term for the majority and minority difference is formulated as ( 1 113 ) and behaves like a sigmoid function. As with the majority carrier mobility, the last two terms in ( 1 112 ) are not functions of temperature and are only used for fitting high PAGE 182 182 impurity concentration data. It should be noted that no extensive experimental data on the minority carrier mobility as a function of temperature is av ailable, according to Klaassen [ 86 ]. Therefore, the temperature fitting parameters were set such that the minority carrier mobility of the proposed model follows the trend of the Philips minority carrier mobility model. The additional parameters required f or fitting the minority carrier data are listed in Table 6 3. The mobility model in equation ( 1 112 ) is compared to experimental data and the Philips model i n Figure 6 6 and Figure 6 7 The comparison is made against the Philips model since it is well formulated for minority carrier mobility. The trend of the proposed model is in agreement with the Philips mobility model for both electron and hole minority car rier mobilities. In order for the majority and minority mobilities to be continuous functions, ( 1 111 ) and ( 1 112 ) the mobilities for electrons and holes can be written as the following set of equations ( 1 114 ) ( 1 115 ) ( 1 116 ) where w is the dopant ratio that allows for the continuous transition between the majority and minority carrier mobilities. Thus mobility as a function of doping levels has been formulated where e,do p defines the electron mobility and h,dop defines the hole mobility. PAGE 183 183 6.3.3 Carrier Carrier Scattering For radiation effects, the carrier carrier scattering effect becomes very important due to the high amount of electron hole pairs that are generated in the de vice during a particle strike. In order to account for carrier carrier scattering, a modified expression of the Conwell Weisskopf formula proposed by Choo [ 87 ] is used and is of the form ( 1 117 ) where n and p are electron and hole densities in cm 3 The doping dependent mobility and carrier carrier scattering mobility terms a re combined using the Mathiessen formula as ( 1 118 ) where the subscript i stands for e or h This results in a unified term for bulk mobility that is a function of doping levels, electron and hole densities, and temperature. The effect of carrier carrier scattering in ( 1 118 ) is compared against experimental data in Figure 6 8. As previo usly discussed, the Philips model highly overestimates mobility at electron hole levels over 10 17 cm 3 I n contrast, the Dorkel Leturcq model uses a similar approach to carrier carrier scattering as the proposed model. The Dorkel Leturcq model fits well for lower carrier concentrations but at concentrations of more than 10 1 7 cm 3 begins to under predict mobi lity. Another issue is that at high injection levels of more than 510 1 9 cm 3 the Dorkel Leturcq model predicts a negative mobility and thus requires an arbitrary minimum mobility condition to be enforced [ 15 ]. PAGE 184 184 In comparison to experimental data, the pro posed model only slightly overestimates the mobility at lower concentrations. However, the electron hole pair concentration generated by a particle strike is typically very high (more than 10 1 7 cm 3 ) [ 15 ]. For this important region, the proposed model cont inues on the assumption that an increase in electron and hole density results in a decrease in carrier mobility [ 87 ]. Above a carrier concentration of 10 17 cm 3 the proposed model predicts a mobility between the Philips and Dorkel Leturcq models and event ually converges to ~2 cm 2 /Vs at a carrier concentration of 10 22 cm 3 Many complications arise when modeling carrier carrier scattering for the ultra high injections conditions that occur following a particle strike. For example, carrier concentrations b ecome degenerate requiring the use of Fermi Dirac statistics, carrier kinetic energies increase, and ambipolar diffusivity increases [ 49 ] Some work has theorized that because carriers are moving together due to ambipolar transport, carrier carrier scatter ing may be minimized suggesting that classical scattering models may not apply for high injection situations [ 88 ] Also, thermalization in the lattice and bandgap narrowing can be factors [ 88 ] Due to these and other complexities, the experimental data sho wn in Figure 6 8 serves as a reminder that more data are needed for carrier concentrations above 10 18 cm 3 6.3.4 Simulation Results and Discussion A series of three dimensional single event transient simulations were run to compare the results obtained using th e proposed mobility model to those obtained from the Philips and Dorkel Leturcq models. The first set of simulation results was also compared to experiments performed by Park et al. [ 42 ] The three mobility models compared in the simulations are the Philip s model, the Dorkel Leturcq model, and the PAGE 185 185 proposed model. A minimum mobility condition (2 cm 2 /Vs) is applied to the Dorkel Leturcq model to prevent the mobility from going negative, as previously discussed. In addition to these three models, a constant m obility model ( e =1417 h = 470.5 cm 2 /Vs) is used to show what occurs when only phonon scattering is consider ed Three different sets of simulations were run to compare the mobility models. In the first set, the mobility models were compared for a 13.5 pJ laser induced current transient and are co mpared to the experimental results that are discussed in detail in chapter 5 [ 42 ] Since the experiment only reached injection levels of 9.810 17 cm 3 two additional sets of simulations were performed to provide insight into the effects of higher injection levels. For the second and third simulation sets, the carrier generation was modeled using a cylindrically symmetrical Gaussian profile more similar an ion strike track. The second set uses the same N + /P diode s tructure as the experiment. For third set, an epitaxial (EPI) N + /P + diode structure was simulated. The simulation variations are summarized in Table 6 4 The dimensions of width, length and depth for the simulation structures were 303040 m and were larg e enough to minimize reflection at the boundaries (Figure 6 9). For each simulation, th e velocity saturation model in equation ( 1 118 ) Shockley Read Hall re combination and Auger band to band recombination models were used. 6.3.5 Experiment Setup The experiment setup for the N+/P diode study is discussed in great detail in chapter 5 so only a brief summary is given here. The diode structure consisted of a heavily d oped n+ region (10 20 cm 3) in a p well (10 18 cm 3 ) that resolved into a p type substrate (10 16 cm 3 ). The n+ and p well junction depths were 0.1 m and 1.5 m, respectively, and a 5 V reverse bias was applied to the device. In the experiment, a PAGE 186 186 cavity dump ed dye laser with a wavelength of 590 nm and a pulse width of 1 ps was used to generate electron hole pairs in the diode (Figure 6 9). The laser direction was normally incident to the diode surface, had a spot size of 12 m in diameter and the energy reach ing the active area of the diode was 13.5 pJ [ 42 ]. 6.3.6 Generated Carrier Distribution For the first simulation set, t he number and distribution of N electron hole pairs generated by the laser pulse was calculated by using the single photon absorption (SPA) equ ation developed by McMorrow as discussed in Chapter 2 [ 22 ] For the second and third simulation sets, the generated electron hole pairs were modeled using a cylindrically symmetrical Gaussian profile. The Gaussian profile had a 1/ e radius of 50 nm, termina ted at a depth of 30 m, and had a linear energy transfer (LET) value of 20 MeV cm 2 /mg. Figure 6 1 0 sh ows the carrier distri bution for the SPA model discussed in Chapter 2 and the cylindrical Gaussian profile. The maximum carrier concentrations for the SPA and Gaussian profiles were 9.810 17 and 1.6410 20 respectively. 6.3.7 Simulation Set 1 Results Experimental Comparison The results of the N+/P diode single event simulations for a laser energy of 13.5 pJ are compared to experimental data in Figure 6 11 and Fi gure 6 12. The simulation result using the proposed model agrees well with the measured data. Data for the experiment were only available up to 10 8 seconds due to the transient measurement setup [ 42 ] As expected, the constant mobility model highly overpr edicts mobility and causes a high current peak and charge collection. The simulation results using the proposed model fall between the Philips and Dorkel Leturcq results. Since the initial maximum electron hole pair concentration is just below 10 18 cm 3 fo r the laser strike, it follows that the proposed model predicts a current transient and charge collection PAGE 187 187 higher than the Dorkel Leturcq model and less than the Philips model due to the high injection mobility shown in Figure 6 8 6.3.8 Simulation Set 2 Results Ion Strike Similar to the previous case, current transients on the N + /P diode due to an ion strike were simulated to provide insight into the effects of higher injection levels. For this set, the cylindrical Gaussian profile in Figure 6 10 was used inst ead of the laser SPA profile. The doping profile and structure are the same as in the previous simulation set. The simulation results of the current transient and charge collection are shown in Figure 6 13 and Figure 6 14. Understandably, the difference in results between the Philips model and the proposed model continues since the difference in high injection mobility increases between the models at higher concentrations (Figure 6 8). The Dorkel Leturcq model still predicts lower charge collection compared to the proposed model. Since the Dorkel Leturcq model underestimates doping dependent mobility (Figure 6 3, Figure 6 6, and Figure 6 7) and predicts lower carrier carrier mobility than the other models (Figure 6 8), it follows that it results in lower cha rge collection than the other models. 6.3.9 Simulation Set 3 Results Epitaxial Diode Current transients for a N + /EPI/P + diode were simulated using the cylindrical Gaussian ion charge deposition profile in Figure 6 10 The diode structure consisted of a heavil y doped n + region (10 20 cm 3 ) on a p type epitaxial substrate (810 14 cm 3 ) placed on a p type substrate (10 20 cm 3 ). The n + junction depth was 0.1 m and the p type EPI layer was 5 m thick. A 5 V reverse bias was applied to the device as in the previous simulations. The simulation results of the current transient and charge collection are shown in Figure 6 15 and Figure 6 16 Due to the much larger depletion PAGE 188 188 region, the charge is collected more quickly than in the case of the bulk diode due to the strong drift region. Once again, the trend continues for charge collection where the simulation results using the proposed model fall between the Philips and Dorkel Leturcq results. 6.3.10 Computational Comparison The proposed model performed well in terms of computati onal efficiency. For example, in a 3 D N+/P diode structure composed of ~6000 volume elements, all device solution times were comparable when separately using each mobility model. The average sum of the matrix assembly and linear solution time was 9.66 sec onds per Newton step for both the Dorkel Leturcq model and the proposed model and 9.73 seconds per Newton step for the Philips model. 6.3.11 Summary A comparison between existing mobility models for device simulation has been presented in section 6.2 to discuss the particular advantages of each model, and a new model (UF high injection mobility model) based on previous formulations is proposed that is computationally efficient and well suited to high injection conditions, such as those found in single event simu lation. As previously discussed, the proposed model has several advantages over the two most commonly used models for radiation effects simulations: the Philips unified mobility model and the Dorkel Leturcq model. The Philips model is formulated in such a way such that it does not match known experimental data for electron and hole concentrations above 10 17 cm 3 The Dorkel Leturcq model was not intended to account for doping concentrations of more than 10 19 cm 3 and was not designed to fit minority mobilit y data. To address the disadvantages of these models, the UF high injection mobility model has been formulated to account for PAGE 189 189 majority and minority carrier mobility, carrier carrier scattering, and temperature dependence making it very suitable for both ra diation effects simulations and general device simulations. Based on the simulation results of both laser and heavy ion charge deposition using the various mobility models the Philips and the Dorkel Leturcq models nt current and charge collection, whereas the proposed model provides an estimate, based on the best data currently available which falls between these bounds. These simulation results indicate that the proposed mobility model gives a peak current, pulse width, and total charge collection for a single event simulation that is closer to experimental measurement than existing mobility models. To aid in mobility model fitting and parameterization, additional experimental data for cases where electron hole car rier densities exceed 10 18 cm 3 will be useful. 6.4 General Purpose Mobility Model Th e previous section describe d a mobility model suitable for single event upset simulations specifically. In this section, a more general purpose mobility will be discussed th at a ccurately describes majority and minority carrier mobilities, carrier carrier scattering, the screening of charge carriers, and temperature dependences. The Philips, previously discussed UF high injection mobility model, and Dorkel Leturcq models acco unt for electron hole scattering in different ways [ 37 ] [ 78], [ 89 ]. It was shown that UF model hold several advantages over the Philips and Dorkel Leturcq models. However, the UF model only focuses on the electron hole scattering mechanism for SEU applicati ons. It does not account for the screening of charge carriers and is dominated by the electron hole scattering component, making it less useful for general purpose device simulation, as shown in the following simulation results section. PAGE 190 190 The focus of the p roposed mobility modeling approach in this section is to accurately fit existing experimental data for lattice, ionized impurity, and electron hole scattering Although some methods formulate mobility starting from fundamental quantum mechanics principles, they can be very computationally intensive and have an adverse effect on simulation time and solution convergence [ 80 ] M obility model s used in device simulation tools start with simplified formulations of lattice and ionized impurity scattering and then use fitting parame ters to match experimental data Our proposed mobility model uses this simplified approach to modeling since finding a balance between physical model accuracy and computational efficiency is important for device simulations. Specifically, the modeling approach in this section uses the Mattheissen rule and follows the same form as equation ( 1 107 ) where bulk silicon mobility can be formulated as ( 1 119 ) with the different components of the mobility represented b y the lattice L donor ND acceptor NA and eh electron hole scattering contributions. In the following subsections, the lattice scattering and majority carrier models are discussed first. Then the minority mobility, electron hole scattering, and char ge screening are defined. Finally, temperature dependence is added to the model and a unified term for mobility is defined. The effect of like carrier scattering (i.e., electron electron, hole hole) is negligible and will be ignored in this study [ 75 ]. PAGE 191 191 6.4.1 La ttice Scattering Carrier scattering in the lattice involves collisions with thermally agitated lattice atoms. The mobility due to this phonon scattering mechanism is a function of temperature and can be written as ( 1 120 ) where the subscript i stands for e (electrons) or h (holes). The mobility dependence on lattice temp erature has been heavily investigated and the parameter is used to fit experimental data [ 82 ]. 6.4.2 Majority Impurity Scattering The majority carrier mobility describes the ionized impurity scattering processes of electrons in n type material (donor sites) a nd holes in p type material (acceptor sites). Our approach to modeling majority mobility is separated into two parts, one for lower doping densities and one for ultra high concentrations. First, the mobility is defined for doping densities below 10 20 cm 3 using the Caughey Thomas model. The Caughey Thomas model is based on plots of experimentally measured mobility data versus the logarithm of doping density, which strongly resemble the Fermi Dirac function [ 84 ]. The Caughey Thomas expression fits experiment al data well for this doping density region and is of the form ( 1 121 ) where C ref and are fitting parameters, N is the doping density, and min and max PAGE 192 192 previously given in ( 1 120 ) as the max term. In an approach similar to Klaassen [ 37 ], the lattice contribution is separated from ( 1 121 ) using the Matthiess en rule and results in the following expression ( 1 122 ) with ( 1 123 ) and ( 1 124 ) where the subscripts ( i, I ) stand for ( e, D ) or ( h, A ) where N D and N A are the donor and acceptor concentrations respectivel y. Experimental data show that mobility drops faster than predicted by the Caughey Thomas expression at concentrations of more than 10 20 cm 3 [ 72 ] This is due to the fact that dopants such as boron, arsenic and phosphorus begin to cluster at higher concen trations [ 82 ] [ 83 ] Since the Caughey Thomas expression no longer matches ( 1 125 ) where C ref and are fitting terms. This formulation uses different fitting terms for electron s and hole s since experimental data show that clustering occurs differ ently depending on dopant type [ 82 ] [ 83 ] The fitting parameters for the majority carri er mobility are given in Table 6 5 To account for the entire range of doping densities, the PAGE 193 193 ionized impurity components in ( 1 122 ) and ( 1 125 ) rule as ( 1 126 ) resulting in a unified term for majority carrier mobility. For example, Figure 6 17 shows electron maj ority mobility in relation to the lattice component in ( 1 120 ) and the ionized impurity components given in ( 1 126 ) An alternate approach to modeling majority carrier mobility would be to use the well known Masetti model formulation [ 72 ]. Although it yields the same results for majority mobility, the Masetti formu lation is not used since the mobility scattering terms in this proposed model are combined strictly by using the Matthiessen rule as in ( 1 108 ) for consistenc y. Since t he Masetti model has been fitted to experimental data very accurately for both electrons and holes, it is compared to our proposed model in Figure 6 18 and Figure 6 19, where it can be seen that the proposed model and the Masetti model agree very well. 6.4.3 Minority Impurity Scattering and Charge Screening Minority carrier mobility is a description of the scattering processes of electrons in p type material and holes in n type material. As with the previous majority carrier formulation, a similar appro ach is used to model the minority carrier mobility by using the Caughey Thomas expression as a starting point. Because minority carrier mobility exceeds majority carrier mobility at high doping concentrations (~110 18 110 20 cm 3 ) [ 86 ] a different set of fitting parameters is used. Following th e modeling approach in ( 1 122 ) the new proposed formulation for minority carrier mobility is of the form PAGE 194 194 ( 1 127 ) where the subscripts ( i, J ) stand for ( e, A ) or ( h, D ) and represents a charge screening parameter. The charge screening parameter is discussed in detail in the next subsection. The mobility model in ( 1 127 ) is compared to both e xperimental data and the Philips model in Figure 6 20 and Figure 6 21 The comparison is made against the Philips model since it is well formulated for minority carrier mobility. The trend of the proposed model is in agreement with the Philips mobility mod el for both electron and hole minority carrier mobilities. The minority mobility fitting parameters for the proposed model are given in Table 6 6. 6.4.4 Electron Hole Scattering and Charge Screening As previously discussed, the electron hole scattering effect be comes very important for radiation effects simulations due to the high density of electron hole pairs that are generated in a device during a particle strike. Similar to the formulation for minority mobility in ( 1 127 ) electron hole scattering is expressed as ( 1 128 ) where the subscripts ( i, K ) stand for ( e, n ) or ( h, p ) where n and p are the electron and hole concentrations respectively and represents the charge screening parameter. The electron hole scattering fitting parameters for the proposed model are given in Table 6 7. The effect of the electron hole scattering in ( 1 128 ) is compared against the Philips model and experimental data in Figure 6 22 As previously discussed, the Philips mobility model is inaccurate at predicting mobility for electron hole densities over 10 17 PAGE 195 195 cm 3 However in comparison to experimental data the proposed model is accurate across the full range of concentrations. For carrier densities below 10 18 cm 3 Dannhauser [ 73 ] and Krausse [ 74 ] measured the sum of electron and hole mobilities as a function of the conce ntration of carriers injected into the weakly doped region of a silicon P I N diode. Unfortunately, for electron hole carrier densities above 10 1 8 cm 3 very few experimental data have been published. However, the proposed model is designed to follow the e xperimental data trend since approximations based on semi classical quantum theory predict that an increase in electron and hole density results in a decrease in carrier mobility [ 15 ]. In terms of radiation effects, the electron hole pair concentration gen erated by a particle strike is typically very high (more than 10 1 8 cm 3 near the center of the particle track ) [ 87 ]. As illustrated by Figure 6 22, it is very important to model this region correctly. A unified mobility term is created by combining the la ttice, majority, minority, and electron expressed as ( 1 129 ) where the subscripts ( i, I, J, K ) stand for ( e, D, A, p ) or ( h, D, A, n ). For example, the electron mobility e,D,A,p is a function of scattering from the lattice e,L donors e,ND accepto rs e,NA and e,p holes. An interesting modeling challenge occurs when establishing expressions for minority and electron hole mobility in a Matthiessen rule scheme as in equation ( 1 120 ) For instance, electron mobility is a function of donor, acceptor, and hole densities as in PAGE 196 196 equation ( 1 129 ) Due to the Matthi essen rule, the mobility term in equation ( 1 129 ) that has the lowest value will dominate the overall mobility value. This behavior becomes a problem for the minority and electron hole components. For example, electron mobility will always be under predicted versus acceptor or hole density since it follows the lowest value for either curve in Figure 6 23 However, the use of the charge screening terms in equations ( 1 127 ) and ( 1 128 ) al lows the mobility to be dominated by the most relevant scattering mechanism. The use of a screening term is valid since at high carrier concentrations carriers tend to screen impurities from other carriers [ 89 ]. The screening terms for electron mobility ar e ( 1 130 ) ( 1 131 ) where the term behaves like a sigmoid function The screening terms indicate that holes screen acceptors just as effectively as acceptors screen holes against electrons. The same assumption is applied to hole mobility, where electrons and donors screen each other. The screening terms for hole mobility are ( 1 132 ) ( 1 133 ) Using the screening term, the electron mobility in p type silicon is determined by the minority mobility term in equation ( 1 127 ) For a particle strike with a high concentration of electron hole pairs, the mobility is dominated by the electron hole PAGE 197 197 scattering term in equation ( 1 128 ) The screening terms allow the proposed model to fit experimental data for both minority mobility and electron hole scattering. Although not physically derived like [ 37 ], the proposed model inherently accounts fo r charge screening in order to fit experimental data. 6.4.5 Temperature Dependence Temperature dependence was previously defined for the lattice scattering limited component of mobility in equation ( 1 120 ) To fit majority carrier mobility to experimental data, two additional fitting terms are added. Rewriting equation ( 1 122 ) as a function of temperature results in the following expression for majority carrier mobility: ( 1 134 ) where T n = T /300 K and the terms are fitting parameters. The temperature fitting parameters are given in Table 6 8. It is important to note that no extensive experimental data on the minority carrier mobility as a function of temperature is available [ 86 ]. Therefore, the temperature fitting parameters were set such that the minority carrier mobility of the proposed model follows the trend of the Philips minority carrier mobility model. A comparison between the prop osed model and measured data for both electron and hole mobilities is given in Figure 6 24 and Figure 6 25 The plots show that the proposed model follows the experimental data trend over a full range of temperatures and doping densities. 6.4.6 Simulation Resul ts Device simulations were run to compare the results obtained using the proposed mobility model to those obtained from other mobility models using the FLOODS PAGE 198 198 simulation tool [ 1 5 ] The three mobility models compared in the simulations are the Philips model the UF model, and the proposed model since they are the most versatile models for general purpose device simulation as shown in Figure 6 1 In addition to these three models, a constant mobility model ( e =1417 h = 470.5 cm 2 /Vs) is used to show what occu rs when only phonon scattering is considered as given by equation ( 1 120 ) For every simulation, the Shockley Read Hall recombination and Auger band to band r ecombination models were used. The simulation results in this work focus on the minority carrier and electron hole scattering components of the mobility models. These are two key areas for the proposed model since accurate experimental data fitting for bot h components is very challenging and has a large impact on simulation results. The minority mobility component is examined in the first set of simulations using a bipolar N/P/N device. For the second simulation set, the electron hole scattering mechanism i s examined using a reverse biased N+/P diode structure. 6.4.6.6 Bipolar N/P/N transistor simulation It is important to model minority carrier mobility accurately for bipolar device simulations. A bipolar N/P/N device serves as a good example since the collector c urrent is due to the injection of electrons from the emitter into the p type base and therefore is a function of the electron minority carrier mobility A set of bipolar device simulations are presented to compare the proposed model versus the Philips mode l. Since minority mobility and charge screening were a focus of the original design, the Philips model provides an excellent and accurate benchmark for a comparison. Additionally, the Philips model was originally designed with bipolar characterization in m ind [ 86 ] Since the focus is to compare mobility models and not to simulate a state of PAGE 199 199 the art device, a very straightforward approach is taken to the N/P/N transistor simulations. The doping profiles of the BJT are represented by step junctions and the di mensions of the device are given in Fig 10. The simulations are performed in 2 D and bandgap narrowing effects are ignored since the focus is mobility modeling. To focus on the minority mobility mechanism for each model, a transient switching simulation is performed. Prior to the transient, the BJT is biased to V BE =1 V and V CE =0.7 V, putting the device into a saturation mode so that the p type base contains a large amount of electron minority carriers (~210 18 cm 3 ). For the transient, V CE remains at 0.7 V and V BE is ramped down from 1 V to 0.3 V (fall time of 1 ps) putting the device into a cut off mode. This voltage switch causes the base to be depleted of electron minority carriers and provides an insightful comparison of how minority mobility modeling affects the device characteristics. As shown by Figure 6 27 the minority mobility component plays a large role in the results. The proposed model agrees well with the Philips model with only a 3% error for the saturation mode current. Since scattering i s minimal for the constant mobility model, the current is highly over predicted when compared against the Philips model. The UF model vastly under predicts current because of the dominant electron hole scattering term which is the focus of the UF model. In saturation mode, the base region of the BJT contains a high number of both electrons and holes. In the UF model, the electron hole scattering is modeled using a modified expression of the Conwell Weisskopf formula proposed by Choo [ 87 ] and is of the prop ortionality of ( 1 135 ) PAGE 200 200 where n and p are electron and hole densities in cm 3 As evident by ( 1 135 ) for any condition in which the electron and hole densities are high, the mobility will be very low. Interestingly, with this scattering term neglected, the UF model was accurate to within 5% of the saturation current predicted by the Philips model. Therefore, although suitable for majority carrier devices such as MOSFETs and for single event simulations where high densities of e lectron hole pairs are prevalent, the UF model is poorly suited for characterizing minority carrier devices such as BJTs. The proposed model does not suffer from this effect due to the formulation of mobility in ( 1 129 ) and electron hole scattering in ( 1 128 ) 6.4.6.7 N+/P diode simulation In this simulation set, t he mobil ity models are compared for a laser induced current transient and are co mpared to experimental results The influence of electron hole since a large number of electron hole p airs are generated along the laser strike path [Cu10]. The experimental and simulation setup will only be briefly described since very detailed descriptions of the experiment and simulation setup are given in chapter 5. In the experiment, a cavity dumped dye laser with a wavelength of 590 nm and a pulse width of 1 ps was used to generate electron hole pairs in the diode (Figure 6 28 ). T he number and distribution of N electron hole pairs generated by the laser pulse was calculated by using the single photon absorption (SPA) equation developed by McMorrow [ 22 ] Th e maximum carrier concentration for the SPA profile w as 9.810 17 cm 3 The results of the N+/P diode single event simulations for laser energy of 13.5 pJ are compared to experimental data in Figure 6 29 and Figure 6 30 Data for the experiment were only available up to 10 8 seconds due to the transient measurement PAGE 201 201 setup [ 42 ] The simulation result using the proposed model agr ees well with the measured data and the UF model. Because the UF model was de signed specifically for SEU simulations, the result shows that the proposed model works for cases of high injection quite well. As expected, the constant mobility model highly over predicts mobility and causes a high current peak and charge collection. Sin ce the initial maximum electron hole pair concentration is just below 10 18 cm 3 for the laser strike, it follows that the proposed model predicts a current transient and charge collection less than the Philips model due to the high injection mobility shown in Figure 6 22 6.4.7 Computational Comparison The proposed model performs well in terms of computational efficiency. For example, in a 3 D N+/P diode structure composed of ~8000 volume elements, all device solution times were comparable when separately using each mobility model. The average sum of the matrix assembly and linear solution time per Newton step was obtained and when compared against the result using the Philips model, the UF model was 3.6% faster and the proposed model was 6.5% faster. 6.4.8 Summary A c omparison between existing mobility models for device simulation was presented in section 6.2 to illustrate the particular advantages of each model, and a new, computationally efficient model based on both previous and new formulations is proposed. The pro posed model is well suited for high injection conditions like those found in SEU simulations and for conditions where minority carrier mobility is important, such as bipolar devices. The proposed model has several advantages over the two most recent models used for radiation effects simulations: the Philips unified mobility model and the UF model. The Philips model is formulated in such a way such that it PAGE 202 202 does not match known experimental data for electron and hole concentrations above 10 1 7 cm 3 Although a ccurate for SEU simulations, the UF model suffers from a dominating electron hole scattering term, making it inaccurate for bipolar transistor simulations. To address the disadvantages of these models, the proposed mobility model has been formulated to fit experimental data for majority and minority carrier mobility, carrier carrier scattering, and temperature dependence. The simulation results show that the proposed model is very suitable for both radiation effects simulations and general purpose device si mulations. 6.5 Interface Mobility Models 6.5.1 Lombardi Model For devices such as MOSFETs, carriers are subjected scattering by acoustic surface phonons and surface roughness at the semiconductor insulator interface. These effects dominate the mobility at the channe l interface, whereas the bulk mobility dominates in low field regions away from the inversion layer. The bulk mobility term in ( 1 129 ) can be used with existi ng models that account for the degradation of mobility at interfaces such as those formulated by Lombardi [ 90 ] and Darwish [ 38 ]. In these approaches, the transverse field E dependent mobility terms are combined with the bulk mobility term using the Matthi essen rule as ( 1 136 ) where b represents the bulk mobility formulate d in ( 1 129 ) ac the acoustic phonon scattering, and sr the surface roughness scattering. Since the interface models are already very well fit to experiment al data, the mobility defined in equation ( 1 136 ) is PAGE 203 203 used as given in [ 38 ]. An example of the mobility dependence on effective field is given in F igures 6 31, 6 32, and 6 33. 6.5.2 Velocity Saturation Model To account for high field saturation, the Canali [ 91 ] approach can be used and is formulated as ( 1 137 ) where 0 is the low field mobility, E  is the driving field, and is a temperature dependent fitting parameter. The Canali model also is based on the Caughey Thomas formula a s in equation ( 1 121 ) and is commonly used in device simulation programs. A plot of electron and hole drift velocity versus electric field is given in Figure 6 34. PAGE 204 204 Table 6 1. Majority Carrier Mobility Fitting Parameters at 300 K (High injection Model) Parameter Electrons (in n type Si) Holes (in p type Si) max 1417 470.5 0 52.2 44.9 1 39.4 29.0 1 0.68 0.719 2 2.0 2.0 C ref,1 9.6810 16 2.2310 17 C ref,2 3.4310 20 6.1010 20 Table 6 2 Temperature Dependence Fitting Parameters (High injection Model). Parameter Electrons Holes 0 0.57 0.57 1 2. 33 2.33 2 2.4 2.4 3 0.4 0.4 4 2.33 2.8 Table 6 3. M inority Carrier Mobility Fitting Parameters (High injection Model). Parameter Electrons (in p type Si) Holes (in n type Si) 2 1270 370 3 39 33 4 150 100 C ref,3 4.6810 16 110 17 C re f,4 3.3410 20 3.3410 20 C ref,5 210 20 210 20 4 3.7 3.7 PAGE 205 205 Table 6 4. Overview of Simulation Variables (High injection Model). Simulation Set Set 1 Set 2 Set 3 Comparison to experiment data Yes No No Structure type N + /P diode N + /P diode Epit axial N + /P Generated electron hole pair profile Single Photon Absorption Energy=13.5 pJ Gaussian LET = 20 MeV cm 2 /mg Gaussian LET = 20 MeV cm 2 /mg Table 6 5 Majority Carrier Mobility Fitting Parameters (General Purpose Model). Parameter Electrons (in n type Si) Holes (in p type Si) max 1417.0 470.5 min 68.5 44.9 1 72.0 49.6 2 1489.0 520.1 3 10 19 4 1417.0 470.5 C ref,1 9.6810 16 2.2310 17 C ref,2 910 19 1.510 20 1 0.711 0.719 2 2 2 Table 6 6 M inority Carrier Mobility Fitting Parameters (General Purpose Model). Param eter Electrons (in n type Si) Holes (in p type Si) 5 525.4 552.7 C ref,3 1.810 17 4.010 17 3 0.6 0.75 0.55 0.55 Table 6 7 Temperature Fitting Parameters (General Purpose Model). Parameter Electrons Holes 2.27 2.25 2 0.1 0.5 3 0.2 0.1 PAGE 206 206 Figure 6 1 Qualitative comparison of commonly used bulk silicon mobility models for device simulation Figure 6 2 Majority electron mobility as a fu nction donor concentration for different mobility models at 300 K. PAGE 207 207 Figure 6 3 Majority hole mobility as a function of acceptor concentration for different mobility models at 300 K. Figure 6 4 Majority electron mobility as a function of temperature and donor concentration. Symbols represent experimental data from [ 82 ]. PAGE 208 208 Figure 6 5 Majority hole mobili ty as a function of temperature and acceptor. Symbols represent experimental data from [ 83 ]. Figure 6 6 Minority electron mobility in p type silicon at 300 K. Symbols represent experimental data fr om Swirhun [ 93 ], Dziewior [ 94 ], Tang [ 95 ]. PAGE 209 209 Figure 6 7 Minority hole mobility in n type silicon at 300 K. Symbols represent experimental data from Dziewior [ 94 ], Burk [ 96 ], Swirhun [ 97 ], Wang [ 98 ]. Figure 6 8 Sum of electron and hole mobility as a function of carrier concentration versus experimental data at 300 K. Symbols represent experimental data from [ 73 ] [ 74 ]. PAGE 210 210 Figure 6 9 Schematic of laser induced current transients [ 42 ] and 3 dimensional simulation structure of the N+/P diode, 303040 m. Figure 6 10 Electron hole pair distributions used in the simulations. A ) Single photon absorption, laser energy = 13.5 pJ, radius = 6 m [ 22 ] B ) Cylindrical Gaussian, LET = 20 MeV cm 2 /mg, 1/ e radius = 50 nm PAGE 211 211 Figure 6 11 Simulated laser induced current transients in a reverse biased Si N+/P diode. Compared to experimental data for a laser energy of 13.5 pJ. Figure 6 12 FLOODS predicted charge collection in a reverse biased Si N+/P diode. Compared to experimental data for a laser energy of 13.5 pJ. PAGE 212 212 Figure 6 13 Simulated current transients in a reverse biased Si N+/P diode. Strike track modeled by a cylindrical Gaussian, LET = 20 MeV cm 2 /mg. Fig ure 6 14 FLOODS predicted charge collection for a reverse biased Si N+/P. PAGE 213 213 Figure 6 15 Simulated current transients in a reverse biased Si N+/EPI/P+ diode Strike track modeled by a cylindrical Gaussian, LET = 20 MeV cm 2 /mg. Figure 6 16 FLOODS predicted charge collection for a reverse biased Si N + /EPI/P + diode. PAGE 214 214 Figure 6 17 Contributions to the majority electron mobility as given by equation ( 1 126 ) Figure 6 18 72 ] for majority electron mobility as a function donor concentration at 300 K. PAGE 215 215 Figure 6 19 model [ 72 ] for majority hole mobility as a function acceptor concentration at 300 K. Figure 6 20 Minority electron mobility in p type silicon at 300 K. Symbols represent experimental data from Swirhun [ 9 3 ], Dziewior [ 94 ], Tang [ 95 ]. PAGE 216 216 Figure 6 2 1 Minority hole mobility in n type silicon at 300 K. Symbols represent experimental data from Dziewior [ 94 ], Burk [ 96 ], Swirhun [ 97 ], Wang [ 98 ]. Figure 6 22 Sum of electron and hole mobility as a function of carrier concentration versus experimental data at 300 K. Symbols represent experimental data from [ 73 ] [ 74 ]. PAGE 217 217 Figure 6 23 Comparison of electron mobility as a function of acceptor site and/or hole density. Figure 6 24 Majority electron mobility as a function of temperature and donor concentration. Symbols repres ent experimental data from [ 82 ]. PAGE 218 218 Figure 6 25 Majority hole mobility as a function of temperature and acceptor concentration. Symbols represent experimental data from [Li78]. Figure 6 26 E 0. 1 0. 4 0. 1 N D 10 20 cm 3 N A 10 18 cm 3 N D 10 17 cm 3 N D 10 19 cm 3 B C Not to scale PAGE 219 219 Figure 6 27 FLOODS 2 D simulation results for a saturation to cut off transient. VBE 1.0 V > 0.3 V, VCE =0.7 V. A B Figure 6 28 A ) Schematic of laser indu ced current transients [ 42 ]. B ) Single photon absorption electron hole pair distribution, laser energy = 13.5 pJ, radius = 6 m [ 22 ] PAGE 220 220 Fi gure 6 29 Simulated laser induced current transients in a reverse biased Si N+/P diode. Compared to experimental data for a laser energy of 13.5 pJ [ 42 ]. Figure 6 30 FLOODS predicted charge collection in a reverse biased Si N+/P diode. Compared to experimental data for a laser energy of 13.5 pJ [ 42 ]. PAGE 221 221 Figure 6 31 Enhanced Lombardi electron mobility model (lines) overlaid on the measured mobility data of (po ints) for several doping values [ 38 ]. PAGE 222 222 Figure 6 32. Enhanced Lombardi hole mobility model (lines) overlaid on the measured hole mobility data of (points) for several doping values [ 38 ]. PAGE 223 223 Figure 6 33. Enhanced Lombardi electron mobility model (lines) overlaid on the measured electron mobility data of several researchers at various temperatures [ 38 ]. PAGE 224 224 A B Figure 6 34. Electron (a) and hole (b) drift velocity in silicon as a functio n of electric field at three different temperatures. The points are the experimental data and the continuous line is the best flitting curve obtained with equation (6 34) [ 90 ]. PAGE 225 225 CHAPTER 7 SUMMARY, CONCLUSIONS AND RECOMMENTATIONS FOR FUTURE WORK 7.1 Summary and Conclu sions A wide range of simulation tool enhancements and physical model improvements have been presented in this work. Each chapter gives a general background overview and analytic explanation of each topic to provide a basis for the work. Simulations were p erformed validate the each simulation tool enhancement and to test the accuracy each new physical model. In Chapter 1, a brief overview of single event effects was given, starting with the historical background. Then the radiation environment was discuss ed where the focus was on particle types and radiation sources. An example of how a particle strike can cause a soft was given and it was shown that CMOS device scaling makes microelectr onics more susceptible to single event upset. The simulations tool challenges for single event effects were discussed and a list of possible tool improvements was given. Lastly, the FLOODS/FLOOPS simulation tool that was used for this work was described. I n Chapter 2, detailed descriptions of the physical mechanisms behind single events were given starting with the electron hole pair generation. The physics of carrier ionization and thermalization were described and equations that model particle strike carr ier generation were discussed. The physics behind charge collection mechanisms such as drift, diffusion and funneling were explained and analytic equations for estimating the total charge collection and current transients were given. Next, the effects of d oping, particle energy, mobility, recombination and bandgap narrowing on single event effects were discussed. PAGE 226 226 In Chapter 3, a finite element approach was described which uses the quasi Fermi levels and electrostatic potential as the solution variables. Th is finite element method differs from the conventional finite volume Scharfetter Gummel approach in that it is not restricted to calculating current along the device grid element edges. The Scharfetter Gummel approach works best if the grid is aligned in t he direction of current flow. However, following a particle strike, the carrier movement is isotropic and thus the finite element approach is better suited for this situation. The simulation results show that the finite element approach is faster and more stable for single event simulations. The focus of this Chapter 4 was on finding ways to reduce simulation time, since SEE simulations are very time intensive. The first section described an adaptive gridding scheme which reduces the number grid points (an d thus simulation time) in real time for a single event transient. The second section will discussed a new diffusive boundary scheme that can be used for the non contacted outer edges of a simulation structure. The boundary scheme allows for a smaller devi ce structure to be used for single event simulations which results in simulation time savings. Both the proposed adaptive grid scheme and diffusive boundary sink were simulated and the results for both show an excellent savings in total simulation time. Ch apter 5 discussed the impact of strained silicon on single event behavior. Because front end process induced strain is used in modern CMOS devices, it is essential to model the change in mobility due to stress. A brief overview of the physics of strained s ilicon was given and then the concepts of linear elasticity, strain, and stress were described. Next, a piezoresistance mobility model was formulated and equations were derived to make it transformable to any silicon orientation. Practical applications of PAGE 227 227 the piezoresistance model were shown, started with a uniaxially strained silicon N+/P diode. The experimental and simulation results agreed well when using the piezoresistance mobility model. Finally, the impact of process induced strain on single event be havior for modern CMOS was investigated. Process and device simulations were performed which show that modern strained silicon technology has a minimal impact on single event characteristics for 45 nm CMOS devices, fabricated on (001) wafers with a channel orientation of <110>. However, the CMOS results gave insight into possible SEE mitigation approaches by using strained silicon technology. In the last section, it was shown that using STI regions to induce stress can result in a much lower charge collecti on and current transient for NMOS devices. Chapter 6 described two new bulk mobility modeling approaches for single event simulations in silicon. The first model focuses on modeling the high injection condition that occurs in a particle strike region. The goal of the high injection mobility model was to formulate a mobility model suitable for radiation effects simulations that accurately describes majority and minority carrier mobilities, carrier carrier scattering, and temperature dependences. The second model takes a more generalized approach to mobility modeling and is very suitable as a general purpose mobility model for device simulations. The second model does better at estimating bipolar current flow and also accounts for charge carrier screening. Bo th models are compared against experimental results and single event simulations were run for each. 7.2 Recommendations for Future Work There are many challenges that remain for the simulation of single event effects. This section briefly discusses a few areas in single event modeling that would benefit from additional research. PAGE 228 228 7.2.1 Carrier Ge neration with Hydrodynamic Transport For this work, th e carriers (electron hole charge cloud) have been entered into the simulation at thermal equilibrium. There is currently much debate as to how to correctly model a particle strike for TCAD tools. High l evel Monte Carlo tools such as MRED, model the strike path as a simple cylindrical Gaussian distribution with an associated LET [ 46 ]. On the other end of the spectrum, atomistic simulators look at the interaction of an ion through a material on atom by ato m basis while accounting for the Coulombic interactions. Between these approaches lie device simulation tools, where it would be useful to simulate the strike process and have the option of associating non equilibrium temperatures with the generated carrie rs (in a computationally efficient manner). Therefore, the effects of hydrodynamic transport modeling for single event simulation should be investigated. Drift diffusion transport does not inherently account for carrier temperature and over estimates impac t ionization. The physics of carrier ionization, strike. Additionally, effects that impact deep submicron devices, such as velocity overshoot, are not well modeled b y the drift diffusion model. As an example, the electron energy balance equation for the hydrodynamic model can be written as (7 1) where S n is the energy flux and W n the energy density. The equations for drift diffusion cur rent density are straightforward as discussed in Chapter III. However, electron current density for the hydrodynamic model can be written as (7 2) PAGE 229 229 where the first term accounts for variations in potential, electron affinity, and bandgap. The remaining terms account for carrier temperature gradients, effective mass, and carrier density [ 21 ]. Due to the complexity of the hydrodynamic approach, drift diffusion transport is still the standard for si ngle event device simulation Si mulation time for the hydrodynamic model is problematic due the amount of simulation variables. This is important since simulation time increases with the number of solution variables k as a function of ~ k 3 [ 32 ]. The full form of the hydrodynamic model con sists of eight partial differential equations With this many solutions variables, it may be prohibitive (with respect to simulation time) to use this model in 3 D single event simulations in the near future. 7.2.2 Bandgap Narrowing The bandgap narrowing mode ls that are commonly available (i.e. Slotboom, del Alamo) are a function of doping levels and were described in Chapter II. Because they do not account for electron and hole densities, the bandgap narrowing in a particle strike region may not be accurate [ 88 ]. The effect of bandgap narrowing in the strike region as a function of carrier densities should be investigated for single event effects. A model exists that formulates bandgap narrowing as a function of doping and carrier densities [ 92 ]. The downside of the model is that in order for it to work in a device simulation tool, only the doping density terms can be used. However, using another approach, there may be a numerically efficient way to account for the electron hole pairs densities. 7.2.3 3 D Adaptive Gridding In Chapter 4, an adaptive gridding scheme was demonstrated. However, at the time of this work, the FLOODS simulation tool is not capable of refining regions in 3 D. It PAGE 230 230 would be interesting to investigate the benefits of adaptive gridding in 3 D si nce simulation times are so much longer. Also, a comparison of adaptive gridding for various discretization methods would be useful in 3 D where it would be expected (based on data in Chapter 3) that the finite element quasi Fermi approach would yield the best results. 7.2.4 Single Event Experiments As stated in Chapter 6, it would be useful to have more data for high injection carrier mobility. Currently, data on carrier mobility is limited to n = p =10 18 cm 3 in literature. If more data could be experimentally ob tained, the mobility models in Chapter 6 could be fit to match the electron hole scattering data. This in turn would result in a higher level of accuracy for single event simulations since mobility is a key factor in results. Experiment results for a n uni axially strained N+/P diode are shown in Chapter 5. Additional experiments for strained Si CMOS devices would be especially useful to compare against the simulation results in Chapter 5. The expectation is that a uniaxially strained <110> MOSFET will show a similar trend to the diode experiment results. For a process induced strained CMOS device, it is expected that the change in collected charge and current would be low, since most of the stress is located at the surface of the device. PAGE 231 231 APPENDIX DERIVATIO N OF TRANSFORMABLE P IEZOCOEFFICIENTS This section goes through the full derivation of a fully transformable (orientation) piezoresistance coefficient matrix. Several references were used as starting points for this derivation [ 55 ], [ 83 ]. However, there i s little published literature on piezoresistance transformation for the entire 6x6 tensor matrix, applicable for all orientations. The unprimed and primed coefficients are shown by the following where the direction old Y axis. Figure A 1. Directional cosine angles The direction cosines are given by the following ( A 1 ) The primed piezoresistance coefficients will be derived using the directional cosines as ( A 2 ) The primed piezoresistance coefficient matrix is of the form PAGE 232 232 ( A 3 ) a The first set to be derived is the ( A 4 ) Summing the twenty one terms gives the following : ( A 5 ) This equation can be further simplified using the following identities: ( A 6 ) ( A 7 ) Subs tituting for the first term in parenthesis for equation 6 1 gives: ( A 8 ) ( A 9 ) Then substituting for the second term in parenthesis for equation 6 1 gives: (A 10 ) Leading to the general iiii PAGE 233 233 ( A 11 ) ( A 12 ) ( A 13 ) ( A 14 ) The second set to be derived is the ( A 15 ) Summing the twenty one terms g ives the following : (A 16 ) This equation can be further simplified using the following identities : (A 17 ) (A 18 ) where substitution yields the following: (A 19 ) Leading to the generalized equation for iij j (A 20 ) (A 21 ) PAGE 234 234 (A 22 ) (A 23 ) The last set to be derived is the (A 24 ) Summing the twenty one terms gives the foll owing : (A 25 ) Where substitution (using previously shown identities) yields the following: (A 26 ) Leading to the generalized equation for ijij (A 27 ) (A 28 ) (A 29 ) (A 30 ) Now a complete set of orientation dependent piezoresistance tensors has been derived. As a sanity check and using Kanda [ 61 ] as a reference, for as in PAGE 235 235 Figure A 1, (a common channel orientation for modern CMOS devices is <110>) the full set of piezoresistance tensors can now be written as: (A 31 ) (A 3 2 ) (A 3 3 ) (A 3 4 ) (A 3 5 ) (A 3 6 ) (A 3 7 ) (A 3 8 ) (A 3 9 ) PAGE 236 236 LIST OF REFERENCES [1] P. Dodd, L. W. Masse ngill, "Basic mechanisms and modeling of single event upset in digital microelectronics," IEEE Trans. on Nuclear Science, vol.50, pp. 583 602, 2003. [2] microelectronics," Materials Today, V ol. 9, No. 6, 2006. [3] induced error 2054, 1982. [4] IEEE Trans. Nucl Sci., vol. 43, pp. 561 575, 1996. [5] NSREC Short Course, 2005. [6] vol. 22, pp. 2675 2680, 1975. [7] particle induced soft errors in dynamic 9, 1979. [8] E. Peterson, "Single Event Analysis and Prediction," IEEE NSREC Short Course, 199 7. [9] Course, 2006. [10] 2006. [11] G. Moore, "Progress in Digital Integrated Electronics," IEEE IEDM Tech Digest, pp.11 13, 1975. [12] R. Ronen, A. Mendelson, K. Lai, L. Shih Lien, F. Pollack, and J. P. Shen, 89, pp. 325, 2001. PAGE 237 237 [13] J. Black, D.B. II, W. Robinson, and DM, "Character izing SRAM Single Event Upset in Terms of Single and Multiple Node Charge Collection," IEEE Trans. Nucl. Sci., vol. 55, pp. 2943 2947, 2008. [1 4 ] P. Dodd, F. Sexton, and P. Winokur, "Three dimensional simulation of charge collection and multiple bit upset in Si devices," IEEE Trans. Nucl. Sci., vol. 41, pp. 2005 2017, 1994. [15] M. E. Law, FLOODS/FLOOPS Manual, University of Florida, 2010. Available: http://www.flooxs.ece.ufl.edu / [16] F. Sexton, "Measuremen t of Single Event Phenomena In Devices and ICs," IEEE NSREC Short Course, 1992. [17] L. Massengill, "Single Event Modeling and Prediction Techniques," IEEE NSREC Short Course, 1993. [18] T. Weatherford, "From carrier to contacts, a review of the SEE cha rge collection processes in devices," IEEE NSREC Short Course, 2002. [19] R. A. Reed, R. A. Weller, R. D. Schrimpf, M. H. Mendenhall, K. M. Warren, L. W. Massengill,"Implications of Nuclear Reactions for Single Event Effects Test Methods and Analysis," I EEE Trans. Nucl. Sci., vol. 53, pp. 3356 3362, 2006. [20] J. F. Ziegler, J. P. Biersack, U. Littmark, The Stopping and Rango of Ion in Matter (SRIM). Pergamon Press, New York, 2009. [21] Synopsys, Sentaurus Device Manual, Version 2007. [22] D. McMorrow, W. T. Lotshaw, J. S. Melinger, S. Buchner, and R. L. Pease, "Subbandgap laser induced single event effects: Carrier generation via two photon absorption," IEEE Trans. Nucl. Sci., vol. 49, pp. 3002 3008, 2002. [23] IEEE Trans. Nucl. Sci., vol. 29, pp. 2024 2031, 1982. [24] for Heavy Ions Incident on n and p pp. 4493 4500, 1983. [25] funneling effect on the PAGE 238 238 collection of alpha particle Electron Device Lett. vol. 2, pp. 103 105, 1981. [26] R. Pierret, Semiconductor Device Fundamentals. Addison Wesley, 1996. [27] J. Slotboom and H. Degraaff, "Measurements of bandgap narrowing in Si bipolar transistors," Solid State Electronics, vol. 19, pp. 857 862, 1976. [28] state mino rity carrier transport parameters in heavily doped n Devices, vol. 34, pp. 1580, 1987. [29] D. Klaassen, J. Slotboom, and H. Degraaff, "Unified apparent bandgap narrowing in n and p type silicon," Solid State Electroni cs, vol. 35, pp. 125 129, 1992. [30] Dependence and Related Features of Radiation Ionization 2038, 1968. [31] ison of Discretization 122, September 2009. [32] Aided Design, vol. 4, pp. 462 471, 1985. [33] signal analysis of a silicon Read 77, 1969. [34] Element Approach to Device 1092, 1983. [35] Comput & Visual Sci., vol. 3, pp. 177 183, 2001. [36] M. R. Pinto, "Comprehensive Semiconductor Device Simulation for VLSI", Ph.D. Thesis, Stanford University, August 1990. [37] I. Model State Electronics, vol. 35, PAGE 239 239 pp. 953 959, 1992. [38] M. Darwish, J. Lentz, M. Pinto, P. Zeitzoff, T. Krutsick, and H. Vuong, "An improved electron and hole mobility model for general purpose device simulation," IEEE Trans. on Electron Devices, vol. 44, p. 1529 1538, 1997. [39] The International Technology Roa dmap for Semiconductors, Modeling and Simulation, 2007. Available: http://www.itrs.net / [40] Scheme for Single s. on Nucl. Sci., Dec 2010. [41] R. E. Bank, A. Weiser, "Some A Posteriori Error Estimators for Elliptic Partial Differential Equations." Mathematics of Computation, vol. 44, pp. 283 301, 1985. [42] H. Park, D. J. Cummings, R. Arora, J. A. Pellish, R. A. Reed, R. D. Schrimpf, D. McMorrow, S. E. Armstrong, U. Roh, T. Nishida, M. E. Law, S. E. Induced Current Transients in Strained Trans. Nuclear Sci., vol. 56, pp. 3203 3209, December 2009. [43] D. A. Antoniadis, I. J. Djo Tempered Bulk http://www mtl.mit.edu/researchgroups/Well / [44] 0.1 mm: How High Will 218, 1997. [45] [46] R. D. Schrimpf, R. A. Weller, M. H. Mendenhall, R. A. Reed, L. W. Massengill, Nuclear Instruments and Methods in Physics Research Section B, vol. 261, pp. 1133 1136, 2007. [47] S. E. Thompson, G. Y. Sun, Y. S. Choi, and T. Nishida, "Uniaxial process induced strained Si: Extending the CMOS roadmap," IEEE Trans. Electron Dev., vol. 53, pp. 1010 1020, 2006. [48] K. Cheng, et. al, "A Highly Scaled High Performance 45nm Bulk Logic PAGE 240 240 CMOS Technology with 0.242 m2 SR AM Cell," IEEE Internation Electron Devices Meeting, pp. 243 246, 2007. [49] S. Sze, Physics of Semiconductor Devices. 3rd Ed., Wiley Interscience, 2007. [50] Y. Sun, S. Thompson, and T. Nishida, Strain Effect in Semiconductors. Springer, 2010. [51] N. Mohta and S. E. Thompson, "Mobility Enhancement: The Next Vector to Extend Moore's Law," IEEE Circuits & Devices Magazine, pp. 18 23, 2005. [52] C. Zhi Yuan, M.T. Currie, C.W. Leitz, G. Taraschi, E.A. Fitzgerald, J.L. Hoyt, mobility enhancement in strained Si n MOSFETs fabricated on SiGe on Device Lett., vol. 22, pp. 321 323, 2001. [53] H. Randall, "Application of Stress from Boron Doping and Other Challenges in Silicon Technology ," Masters Thesis, University of Florida, 2005 [54] D. Logan, Finite Element Method. 4th Ed., Thompson, 2007. [55] R. E. Newnham. Properties of Materials. Oxford Univ. Press, 2005. [56] S. E. Thompson, G. Sun, K. Wu, J. Lim, and T. Nishida, "Key differ ences for process induced uniaxial vs. substrate induced biaxial stressed Si and Ge channel MOSFETs," IEDM Tech. Dig., pp. 221 224, 2004. [57] O. A. Amusan, "Analysis of single event vulnerabilities in a 130 nm CMOS technology," Masters Thesis, Vanderbilt University, 2006. [58] J. S. Melinger, S. Buchner, D. McMorrow, W. J. Stapor, T. R. Weatherford, and A. B. Campbell, "Critical Evaluation of the Pulsed Laser Method for Single Event Effects Testing and Fundamental Studies," IEEE Trans. Nucl. Sci., vol. 4 1, pp. 2574 2584, 1994. [59] M. Amiotti, A. Borghesi, G. Guizzetti and F. Nava, "Optical Properties of Polycrystalline Nickel Silicides," Phys. Rev. B, vol. 42, pp. 8939 8946, 1990. [60] C. S. Smith, "Piezoresistance Effect in Germanium and Silicon," Phy s. Rev. vol. 94, pp. 42 49, 1954. PAGE 241 241 [61] Y. Kanda, "A Graphical Representation of the Piezoresistance Coefficients in Silicon," IEEE Trans. Electron Dev., vol. 29, pp. 64 70, 1982. [62] S. E. Thompson et al., "A 90 nm logic technology featuring 50 nm str ained silicon channel transistors, 7 layers of Cu interconnects, low k ILD, and 1um 2 SRAM cell," IEDM Tech. Dig., pp. 61 64, 2002. [63] V. Chan et al., "High speed 45 nm gate length CMOSFETs integrated into a 90 nm bulk technology incorporating strain e ngineering," IEDM Tech. Dig., pp. 3.8.1 3.8.4, 2003. [64] J. Lim, S. Thompson, and J. Fossum, "Comparison of threshold voltage hifts for uniaxial and biaxial tensile stressed n MOSFETs," IEEE Electron Dev. Lett., vol. 25, pp. 731 733, Nov. 2004. [65] K. Mistry et al., "A 45nm Logic Technology with High k+Metal Gate Transistors, Strained Silicon, 9 Cu Interconnect Layers, 193nm Dry Patterning, and 100% Pb free Packaging," IEEE International Electron Devices Meeting, pp. 247 250, 2007. [66] T. Miyashita et al., "High Performance and Low Power Bulk Logic Platform Utilizing FET Specific Multiple Stressors with Highly Enhanced Strain and Full Porous Low k Interconnects for 45 nm CMOS Technology," IEEE International Electron Devices Meeting, vol. 6, pp. 251 254 2007. [67] N. Shah, "Stress Modeling of Nanascale MOSFET," Masters Thesis, University of Florida, 2005. [68] K. Su, Y. Sheu, C. Lin, S. Yang, Wen Jya Liang, Xuemei Xi and A.C. Jaw Kang Her, Yu Tai Chia, Carlos H. Diaz, "A scaleable model for STI mechan ical stress effect on layout dependence of MOS electrical characteristics," Proc. of the IEEE 2003 Custom Integrated Circuits Conference, pp. 245 248, 2003 [69] M. Miyamoto, H. Ohta, Y. Kumagai, Y. Sonobe, K. Ishibashi, and Y. Tainaka, "Impact of reducing STI induced stress on layout dependence of MOSFET characteristics," IEEE Trans. on Electron Devices, vol. 51, pp. 440 443, 2004. [70] R. Arghavani, Z. Yuan, N. Ingle, K. Jung, M. Seamons, S. Venkataraman, V. Banthia, K. Lilja, P. Leon, G. Karunasiri, S. Yoon and A. Mascarenhas, "Stress management in sub 90 nm transistor architecture," IEEE Trans. Electron Devices, vol. 51, p. 1740 1744, 2004. PAGE 242 242 [71] Y. Luo and D. Nayak, "Enhancement of CMOS Performance by Process Induced Stress," IEEE Trans. on Semiconduct or Manufacturing, vol. 18, pp. 63 68, 2005. [72] G. Masetti, M. Severi, and S. Solmi, "Modeling of carrier mobility against carrier concentration in arsenic phosphorus and boron doped silicon," IEEE Trans. Electron Devices, vol. 30, pp. 764 769, 1983. [73] F. Dannhauser, "Dependence of carrier mobility in silicon on the concentration of free charge carriers I," Solid State Electron., vol. 15, p. 1371 1375, 1972. [74] J. Krausse, "Dependence of carrier mobility in silicon on the concentration of free charge carriers II," Solid State Electron., vol. 15, pp. 1377 1381, 1972 [75] B. K. Ridley, Quantum Process in Semiconductors. 2nd Ed., Clarendon Press, pp. 141 152, 1988. [76] T. Zhang, "Impact of Charge Collection Mechanisms in Single Event Effects in SiGe HBT Circuits and Hardening Implications," M.S. Thesis, Auburn University, 2009. [77] A. Balasubramanian, "Measurement and Analysis of Single Event Induced Crosstalk in Nanoscale CMOS Technologies," Ph.D. Dissertation, Vanderbilt University, 2008. [7 8] J. Dorkel and P. Leturcq, "Carrier mobilities in silicon semi empirically related to temperature, doping and injection level," Solid State Electron., vol. 24, pp. 821 825, 1981. [79] P. Chapman, O. Tufte, J. Zook, and D. Long, "Electrical properties of heavily doped silicon," Journal of Applied Physics, vol. 34, pp. 3291 3295, 1963. [80] M. Fischetti, "Effect of the electron plasmon interaction on the electron mobility in silicon," Physical Review B, vol. 44, 1991, p. 5527 5534. [81] N. Arora, J. Haus er, and D. Roulston, "Electron and hole mobilities in silicon as a function of concentration and temperature," IEEE Trans. Electron Devices, vol. 29, pp. 292 295, 1982. [82] S. Li and W. Thurber, "The dopant density and temperature dependence of electron mobility and resistivity in n type silicon," Solid State Electron., vol. 20, PAGE 243 243 pp. 609 616, 1977. [83] S. Li, "The dopant density and temperature dependence of hole mobility and resistivity in boron doped silicon," Solid State Electron., vol. 21, pp. 1109 1 117, 1978. [83] R. F. Tinder, Tensor Properties of Solids. Morgan & Claypool, 2008. [84] D. Caughey and R. Thomas, "Carrier mobilities in silicon empirically related to doping and field," Proc. Inst. Elec. Eng., vol. 55, pp. 2192 2193, 1967. [86] D. Kl aassen, "Physical Modelling for Bipolar Device Simulation," Conference on Simulation of Semiconductor Devices, vol. 4, pp. 23 43, 1991. [87] C. Choo, "Theory of a forward biased diffused junction PLN rectifier: 2," IEEE Trans. Electron Devices, vol. 19, p p. 954 966,1972. [88] J. Laird, S. Onoda, T. Hirao, L. Edmonds, and T. Ohshima, "The Role of Ion Track Structure on High Injection Carrier Dynamics in High Speed Si and III V Optoelectronic Sensors," IEEE Trans. Nucl. Sci., vol. 54, pp. 2384 2393, 2007. [89] D. Cummings, A. F. Witulski, H. Park, R. D. Schrimpf, S. E. Thompson and M. E. Law, "Mobility Modeling Considerations for Radiation Effects Simulations in Silicon," IEEE Trans. Nucl. Sci., vol. 57, no.4, pp.2318 2326, 2010. [90] C. Lombardi, S. Man zini, A. Saporito, and M. Vanzi, "A physically based mobility model for numerical simulation of nonplanar devices," IEEE Trans. Computer Aided Design of Integrated Circuits and Systems, vol. 7, pp. 1164 1171, 1988. [91] C. Canali, G. Majni, R. Minder, and G. Ottaviani, "Electron and hole drift velocity measurements in silicon and their empirical relation to electric field and temperature," IEEE Trans. Electron Devices, vol. 22, pp. 1045 1047, 1975. [92] temperature full random phase app roximation model of vol. 84, pp. 3684 3695, 1998. [93] S. Swirhun, Y. Kwark, R. Swanson, "Measurement of electron lifetime, electron mobility and band gap narrowing in heavily doped p type silicon," PAGE 244 244 IEDM, pp. 2 5, 1986. [94] J. Dziewior and D. Silber, "Minority carrier diffusion coefficients in highly doped silicon," Appl. Phys. Lett., vol. 35, pp. 170, 1979. [95] D. Tang, F. Fang, M. Scheuer mann, T. Chen, and G. Sai Halasz, "Minority carrier transport in silicon," IEDM, pp. 20 23, 1986. [96] D. Burk and V. De La Torre, "An empirical fit to minority hole mobilities," IEEE Electron Device Letters, vol. 5, pp. 231 233, 1984. [97] S. Swirhun, J .D. Alamo, and R. Swanson, "Measurement of hole mobility in heavily doped n type silicon," IEEE Electron Device Letters, pp. 168 171, 1986. [98] C. Wang, K. Misiakos, and A. Neugroschel, "Minority carrier transport parameters in n type silicon," IEEE Tran s. Electron Devices, vol. 31, pp. 1314 1322, 1990. PAGE 245 245 BIOGRAPHICAL SKETCH Daniel Joseph Cummings was born in 1983 in Los Angeles, California. He has three younger siblings and loving parents, Lois and Joseph. He received his B.S. and M.S. in electrical engi neering at the University of Florida in 2006 and 2008 respectively. His hobbies include playing various instruments, cooking, gardening, hiking, and watching college football. Upon completion of his graduate work, he will start working for Intel in Austin Texas. 