Simulat ion Studies of C a r b o n Nanotube Field-Effect Transistors by David Llewellyn John B.A.Sc, The University of British Columbia, 2002 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF Doctor of Philosophy in The Faculty of Graduate Studies (Electrical and Computer Engineering) The University Of British Columbia July 2006 © David Llewellyn John 2006 A b s t r a c t Simulation studies of carbon nanotube field-effect transistors (CNFETs) are presented using models of increasing rigour and versatility that have been systematically developed. Firstly, it is demon-strated how one may compute the standard tight-binding band structure. From this foundation, a self-consistent solution for computing the equilibrium energy band diagram of devices with Schottky-barrier source and drain contacts is developed. While this does provide insight into the. likely be-haviour of CNFETs, a non-equilibrium model is required in order to predict the current-voltage relation. To this end, the effective-mass approximation is utilized, where a parabolic fit to the band structure is used in order to develop a Schrodinger-Poisson solver. This model is employed to predict both DC behaviour and switching times for CNFETs, and was one of the first models that captured quantum effects, such as tunneling and resonance, in these devices. In addition, this model has been used in order to validate compact models that incorporated tunneling via the WKB approximation. A modified WKB derivation is provided in order to account for the non-zero reflection of carriers above a potential energy step. In order to allow for greater flexibility in the CNFET geometries, and to lift the effective-mass approximation, a non-equilibrium Green's function method is finally developed, which uses an atomistic tight-binding Hamiltonian to model doped-contact, as opposed to Schottky-barrier-contact, devices. This approach benefits by being able to account for both inter-and intra-band tunneling, and by utilizing a quadratic matrix equation in order to improve the computation time for the required self-energy matrices. Within this technique, an expression for the local inter-atomic current is derived in order to provide more detailed information than the usual compact expression for the terminal current. With this final model, an investigation is presented into the effects of geometrical variations, contact thicknesses, and azimuthal variation in the surface potential of the nanotube. Table of Contents Abstract 1 1 Table of Contents i U List of Figures v l Acknowledgements v m Dedication l x 1 Introduction 1 1.1 Background 1 1.2 Simulation Studies 4 2 Band Structure 6 2.1 Graphene 6 2.2 Zone-Folding 9 3 Equilibrium 13 3.1 Electrostatics 1 4 3.1.1 Closed Geometry 14 3.1.2 Open Geometry , 15 3.2 Solution Procedure 16 3.3 Results and Discussion 17 4 The Effective-Mass Method 21 4.1 Finite Difference Electrostatics 22 iii 4.2 Treating Quantum Phenomena 23 4.2.1 Asymptotic Treatment of Tunneling 23 4.2.2 Solving Schrodinger's Equation 26 4.2.3 Quantum Capacitance 30 4.3 Planar Devices • 34 4.3.1 Finite Element Electrostatics 35 4.4 Solution Procedure 37 4.5 Results and Discussion 39 4.5.1 D C 39 4.5.2 Switching Speed 43 5 T h e T i g h t - B i n d i n g M e t h o d 51 5.1 Relevance of the Green's Function 51 5.2 The Hamiltonian 55 5.3 Computing the Green's Function 56 5.3.1 Green's Function for a Semi-Infinite Lead 59 5.3.2 Test Cases 6 1 5.4 Current Operators 64 "5.4.1 Effective-Mass 65 5.4.2 Tight-Binding • • . 66 5.4.3 Current Operator Summary 70 5.5 Alternate Current Expression 72 5.6 Solution Procedure 74 5.7 Results and Discussion 78 6 C o n c l u s i o n s 87 6.1 Current Progress • • 87 6.2 Future Work 90 B i b l i o g r a p h y 92 iv A Solution to Laplace's Equation 107 A. l Potential Applied to the Curved Surface 108 A.2 Potential Applied to One End 110 A.3 Total Solution 112 B Green's Function for a Cylindrical Box 113 C Conformal Mapping of the Cylindrical Laplace Equation 117. C.l Transformation of the Cartesian Laplace Equation 118 C. 2 Transformation of the Cylindrical Laplace Equation 120 D The W K B Approximation 122 D. l Approximate Solution • • • • 123 D.2 Error 123 D. 3 Turning Points and Matching 124 E Solving for the Newton Update 129 E. l Computing the Raw Newton Update 129 E.2 Using Exact Line Searches 130 v L i s t o f F i g u r e s 1.1 Schematic graphene lattice 2 2.1 Dispersion relation and density of states for a (16,0) carbon nanotube . . 12 3.1 Schematics of the closed and open cylindrical CNFET geometries 13 3.2 Equilibrium energy bands near the source for a 5.6 nm gate-oxide thickness 19 3.3 Equilibrium energy bands near the source for a 30.9 nm gate-oxide thickness . . . . 19 3.4 Equilibrium energy bands with a permittivity five times higher than in Fig. 3.2 . . . 20 4.1 Model energy band diagram 24 4.2 A parabolic fit to the (16,0) CN band structure 27 4.3 ID quantum capacitance at equilibrium 32 4.4 ID quantum capacitance for applied bias voltages 35 4.5 Planar geometry CNFET 36 4.6 ID piecewise linear basis function 38 4.7 Conduction band edges for VQS = 0.5 V, and VDS = 0 and 0.4 V 40 .4.8 Conduction band edges and transmission probabilities 41 4.9 Carrier concentration as a function of energy 42 4.10 TD-VGS for a closed cylindrical device 42 4.11 Schematic showing the ambipolar conduction mechanism 43 4.12 Circuit for computing the delay time • • • • 44 4.13 Computation of the gate capacitance 45 4.14 Cross-section of the planar geometry CNFET 46 4.15 1D-VGS for three planar devices 47 vi 4.16 Intrinsic switching time, and ON/OFF ratio 48 4.17 CN charge as a function of the gate voltage 49 4.18 Valence bands for two oxide thicknesses 50 4.19 Extrinsic switching time, and ON/OFF ratio 50 5.1 Tight-binding density of states using the self-energy method 61 5.2 Local density of states for longitudinal and azimuthal potential variation 63 5.3 Sample one- and two-dimensional lattices 67 5.4 Hexagonal lattice with a primitive unit cell 69 5.5 Atomic structure of a (2,0) carbon nanotube 72 5.6 Schematics of the double-gate and semi-cylindrical-gate CNFET geometries 78 5.7 Local CN potential for the double-gate and the semi-cylindrical-gate devices . . . . 79 5.8 Local CN potential as a function of azimuthal angle 80 5.9 Local density of states and energy spectrum of the electron density 81 5.10 Energy spectrum of the current 82 5.11 TZ>-VGS for double-gate and semi-cylindrical-gate devices 83 5.12 Energy band diagrams for the semi-cylindrical-gate geometry 84 5.13 Inter-atomic current for a double-gate device 85 5.14 Comparison of NEGF results to those using the Flietner dispersion relation 86 vii A c k n o w l e d g e m e n t s I have been very fortunate to work with an outstanding group over the course of my studies. In particular, I am grateful for early discussions with Paul Pereira and Jason Clifford, and more recently, with Dylan McGuire. I have also benefited from many detailed discussions with Leonardo Castro, as we simultaneously pursued a deeper understanding of carbon nanotube devices. Finally, it has been a privilege to work with my mentor, Professor David Pulfrey, whose guidance, expertise, and insight have been vital to this research. viii D e d i c a t i o n This work would not have been possible without the love and support of my family. My children, Alicia, Emlyn, and Arden, are constant sources of wonder and inspiration. In addition, I owe everything to my wife, Kerrene, who has contributed to this work in more ways than I can list here. She is an inspiration in all of the aspects of my life, and this thesis is as much a result of her hard work and dedication as it is mine. I dedicate this work, and much more, to them. ix Chapter 1 I n t r o d u c t i o n The aggressive scaling of the metal-oxide-semiconductor field-effect transistor (MOSFET) to ever-smaller physical dimensions has required the introduction of a wide variety of technologies in order to fabricate these devices [1]. With each generation, however, the manufacturing challenges are becoming increasingly difficult, and according to the International Technology Roadmap for Semi-conductors [1], intensive research is needed in order to continue this process, and indeed, to develop novel devices and methods that will move the technology improvements in other directions. A candidate transistor that may allow for both the shrinking process to continue, and for the de-velopment of novel architectures, is the carbon nanotube field-effect transistor (CNFET). Essentially, this device is a conventional MOSFET, except that its channel is made up of a carbon nanotube (CN), and that the source and drain interfaces may be Schottky contacts. In the following section, a brief background to CNs is provided, with reference to some of their interesting properties. 1.1 Background In 1976, the growth of cylindrical tubes of carbon was reported by Oberlin et al. [2], but they were not recognized as having the CN structure, and were not studied in detail. Intensive interest in CNs did not begin until 1991 when Iijima noticed their growth on an arc discharge electrode [3], and used electron microscopy to probe their atomic structure. Since that time, a great deal of theoretical and experimental work has gone into understanding carbon nanotubes. To begin, consider the CN atomic structure. Conceptually, it may be thought of as a graphene sheet that has been rolled into a seamless cylinder. Shown in Fig. 1.1 is a schematic of the hexagonal lattice of carbon atoms that form the atomic monolayer of graphene. The two primitive translation vectors, ai and a2, point to equivalent sites in the lattice, i.e., points that differ by integer multiples of these two vectors are identical in an infinite graphene sheet. Considering two nearest-neighbour 1 atoms, primitive translations may be performed in order to fill the entire graphene sheet with no overlaps or gaps. As a result, graphene may be defined by simply specifying the primitive translation vectors and the two-atom basis. In order to form a seamless cylinder, a patch is taken out of the infinite sheet, and rolled such that the circumference is defined by the chiral vector L = ruai + 712^2, where n\ and ni are integers. This maps an atom in the lattice to coincide with an equivalent point specified by L . Figure 1.1: Schematic graphene lattice showing the primitive translation vectors. Carbon nanotubes are uniquely defined, then, by simply specifying the (711,712) indices. Note that, in Cartesian coordinates, (1.1) 3a „ a\/3 „ a i = T x + ~ y ' 3a „ a\[2> „ a 2 = Yx~— y' (1.2) where a is the nearest-neighbour carbon-carbon spacing of about 1.4 A [4], and x and y are unit vectors in the x- and y-directions, respectively. In terms of some basic terminology, CNs are named "armchair" if 711 = Ti2, "zigzag" if 712 = 0, and "chiral" otherwise. The armchair and zigzag names relate to the carbon-carbon bond pattern that may be observed when looking along the circumference of the tube. An interesting property of CNs is that they may be metallic or semiconducting depending on the choice of (711,712). Specifically, the tube is metallic if n\ — 712 is a multiple of 3, and is semiconducting otherwise with a bandgap, EQ, given approximately by [4-6] (1.3) where 7 is an overlap integral which will be defined in Chapter 2, and Rc is the CN radius. Since semiconducting tubes are possible, this raises the prospect of using CNs in order to make semicon-ductor devices, such as transistors. Indeed, this was performed in 1998 by Tans et al. [7], and this early device stimulated a great deal of interest in producing, and simulating, reliable CNFETs. Perhaps one of the most appealing features of CNs, for nanoelectronics, is their near-ballistic transport due to a limited carrier-phonon interaction [8]. Mean-free-paths for acoustic phonons have been reported in the range of 300 nm to greater than 2.9 /zm [9-12], and for optical phonons the range is typically in the 10-100 nm range [12,13]. In Ref. [14], it is asserted that transport should be ballistic in CNs up to 60 nm in length, while Ref. [10] gives a value of 20 nm for a CN in a low field. Mobilities of 104-105 cm2/Vs have been reported [9,13], and a recent theoretical prediction also falls within this range at low field [15]. When compared to the bulk silicon values of around 102 and 103cm2/Vs for holes and electrons [16], respectively, the CNFET appears to hold some promise. For silicon nanowire transistors, the comparison favours CNs even more heavily since the average mobility is in the range 30-560 cm2/Vs, with a maximum value of 1350cm2/Vs [17]. Couple this with the large current densities that CNs can withstand [18,19], roughly three orders of magnitude greater than that reported for silicon nanowires [20], and the carbon nanotube looks even more promising. In addition, with their cylindrical shape, CNs may lend themselves to the realization of a coaxial field-effect transistor, which is predicted to be the ultimate geometry in terms of good short-channel behaviour [21-23]. This structure is being pursued experimentally [24], and has been successfully made using vertical nanowires [25], but not yet with carbon nanotubes. High permittivity dielectrics are also a possibility towards reducing short-channel effects due to the increased capacitative coupling between the channel and the gate, and success has been achieved using Zr02 [26] and SrTi03 [27], with relative permittivities of 25 and 175, respectively. High permittivities have also been achieved using electrolytes [28,29], where the relative permittivity of the solution could be around 80 [30], although liquid electrolytes are likely to have limited utility in highly integrated circuits. Furthermore, there may be biological applications of CN devices [31], for example, in ultrafast DNA sequencing or as a DNA-incorporated switch [32]. Structures are also being investigated for use as chemical sensors [33,34]. Additionally, with the unique electrical properties of CNs, one could envision all-carbon electronics with metallic nanotube leads connected to semiconducting nanotube devices [35,36]. 3 Given these properties, it is clear why there is interest in studying carbon nanotubes and carbon nanotube electronics. There are, however, some unresolved problems in CNFET fabrication. Early devices were fabricated by dispersing CNs from CN ropes, and then positioning them using an atomic force microscope tip [37]. While this method has been largely replaced by the more reliable chemical vapour deposition (CVD) technique [38,39], the major obstacle of chirality control remains an outstanding problem. Chirality control is important as this determines whether the CN is metallic or semiconducting, and if semiconducting, it also determines the energy bandgap. Since chirality control is difficult, another option is to devise methods to separate CNs based on chirality or conductivity. This may be facilitated by the presence of single-stranded DNA [31]. Sorting of CNs by diameter [40], by metals and semiconductors [41], and by chirality [42,43], have been reported, and recent work has shown that semiconducting tubes may be preferentially grown in the plasma enhanced CVD method [44], with a yield of around 90% semiconducting CNs. With this latter technique, perhaps it would be sufficient to control only the diameter, and hence the bandgap. This may prove to be a more attainable goal. 1.2 Simulation Studies In order to examine the expected performance, and performance limits, of CNFETs, detailed sim-ulations of these devices are presented. In Chapter 2, the band structure and density of states are computed. This serves as the foundation for this work. Chapters 3-5 present a range of CNFET models and results: beginning with a simple case, invoking approximations that reduce the problem complexity, and then successively stripping away some of those simplifying assumptions such that each subsequent model captures more of the device physics. In particular, Chapter 3 considers the equilibrium case of devices with Schotkky-barrier contacts at the source and drain [45-50]. This is the simplest case to consider under the assumption that potential energy variations are sufficiently gradual to allow for the usual Fermi-Dirac statistics to apply. As will be seen in Chapter 4, rapid potential variations may lead to resonance effects which serve to alter the distribution function. These equilibrium simulations resulted in the development of a compact model for CNFETs that was an early attempt at including the effect of tunneling through the Schottky barriers [51]. For steep potential energy variations, and for the non-equilibrium case, Chapter 4 presents a self-4 consistent Schrodinger-Poisson solver, which uses the effective-mass approximation to the CN band structure. At the time of this work, CNFETs were being modeled using electrostatics appropriate for an infinite CN, and under the assumption that transport was limited by thermionic emission over the top of a potential energy barrier [26,52]. Here, a modified WKB expression is contributed that is appropriate for use in compact models [53], and that accounts for the reflection of carriers above a potential energy step, such as in a metal-semiconductor contact where the Fermi wavevector changes abruptly upon crossing the interface [54]. It transpired that the new modified WKB expression was very successful in reconciling the compact model results with those of the, more detailed, full effective-mass model [53]. While a cylindrical structure is required in order for this method to be strictly valid, it may also be reasonable for less symmetric structures if the azimuthal potential energy variation is small enough [55]. In order to examine such variation, which had apparently been neglected in order to speed up similar calculations [56], and to lift the effective-mass approximation, Chapter 5 details a non-equilibrium Green's function method for computing the local density of states from Schrodinger's equation. This is, again, coupled with Poisson's equation in a self-consistent fashion, but here an atomistic, nearest-neighbour tight-binding Hamiltonian, as opposed to the previously employed effective-mass variety, is used. This method was developed independently of that due to Fiori et al, despite the similarities in technique [57,58]. One key difference is in the computation of the self-energy matrices, where a different algorithm is presented herein. In addition, expressions for the inter-atomic current are derived, and these may prove useful when studying CNs with defects. The focus in this chapter is on doped-contact devices without Schottky-barriers at the end contacts, as may be appropriate for some devices [59-62]. The p-i-n case [63-65] is discussed as it may exhibit less than a 60mV/dec subthreshold slope, which is the limit for devices governed by the thermionic emission of carriers over a barrier. A lower subthreshold slope would be a tremendous boost for highly integrated logic circuitry, which is presently curtailed by the static power dissipation from supposedly "OFF" gates. Chapter 6 presents a summary of contributions, key conclusions to be drawn from this work, and outlines future work that may be built from the solid foundation presented herein. In addition, key derivations are included in the appendices in order to facilitate the reproduction, and extension, of this work by other investigators. 5 Chapter 2 Band Structure In order to simulate a CNFET, the band structure of a carbon nanotube must first be determined. To this end, the nearest-neighbour, tight-binding approximation is applied. From the band structure, the density of states (DOS) may be calculated, and this is required in order to compute the amount of charge on the tube. This fundamental calculation is a basic tool, and is an important first step to take in understanding the behaviour of CNFETs. 2.1 Graphene Consider, first, the band structure of graphene. The overall periodicity of the lattice is defined by a two-point basis which may be shifted by integer multiples of the primitive translation vectors in order to tile, with no overlaps or gaps, an infinite graphene crystal. The electron wavefunctions are governed by the time-independent Schrodinger equation, HV{r) = £ # ( r ) , (2.1) where H is the Hamiltonian, 3* is the wavefunction, r is position, and E is energy. Following the usual tight-binding procedure [66], the Hamiltonian may be decomposed into H = Hat + AU(r), where Hat is the Hamiltonian for an isolated atom, and AU(r) represents the modification to this Hamiltonian due to the presence of all of the other atoms. The atomic Hamiltonian has a complete set of eigenfunctions, ^ a t . i i with their associated eigenenergies, Ei, i.e., i ? a t ^ a t , t = ^ t ^ a t . i - Further, due to the underlying periodicity of the graphene lattice, Bloch's theorem [66] may be used to write the total wavefunction as *(r) = 53 e i k R [BiS(r - R) + B2E(r - d - R)], Ft 6 w h e r e t h e s u m m a t i o n i s o v e r a l l r e c i p r o c a l l a t t i c e v e c t o r s R, k i s t h e w a v e v e c t o r , d i s a v e c t o r f r o m o n e a t o m i n t h e t w o - p o i n t b a s i s t o t h e o t h e r , B\ a n d £ ? 2 a r e c o n s t a n t s , a n d S i s k n o w n a s a W a n n i e r f u n c t i o n [66]. S i n c e t h e $ a t , » a r e a c o m p l e t e s e t , S m a y b e w r i t t e n a s S(r) = ^Ci*a t , i ( r ) , i w h e r e Ci a r e c o n s t a n t s . U p u n t i l t h i s p o i n t , t h e e x p r e s s i o n s h a v e b e e n e x a c t . T h e f i r s t a p p r o x i m a t i o n i s t h a t t h e W a n n i e r f u n c t i o n s a r e w e l l - d e s c r i b e d b y a s i n g l e c h o i c e o f e i g e n f u n c t i o n . N a m e l y , l e t Ci = 1 f o r t h e pz o r b i t a l s , d e f i n e d b y i = 1 , a n d C7» = 0 o t h e r w i s e . M u l t i p l y i n g t h r o u g h b y VP^ t -(r) a n d i n t e g r a t i n g a l l o w s E q . ( 2 . 1 ) t o b e w r i t t e n a s J * i t i J . ( r ) t f a t * ( r ) d 3 r + J ^ t i J . ( r ) A £ / ( r ) * ( r ) d 3 r = J *{tJ(r)E<i>(r) d 3 r , ( 2 . 2 ) w h e r e t h e i n t e g r a l i s o v e r a l l s p a c e , ^ a t j i s a n a r b i t r a r y e i g e n f u n c t i o n , a n d f d e n o t e s t h e H e r m i t i a n c o n j u g a t e . N o t e , h e r e , t h a t E i s a s c a l a r , s o i t m a y s i m p l y c o m e o u t s i d e t h e i n t e g r a l . I n a d d i t i o n , s i n c e t h e a t o m i c H a m i l t o n i a n o p e r a t o r i s H e r m i t i a n , j(r)Hat = (Hat^>atj(r))^. S u b s t i t u t i n g i n t h e a b o v e d e f i n i t i o n s y i e l d s (E - Ej) / 5 ~ J e i k - R ^ t J . ( r ) [ B i * a t i i ( r - R ) - l - B 2 « ' a t . i ( r - d - R ) ] d 3 r R = /Ee i k R* a t J(r)At/(r) [Bi* at,i(r - R) + B 2 *at,i(r - d - R)] d 3 r . R C o n s i d e r i n g t h e c a s e w h e r e j = 1 , a n d r e a r r a n g i n g , g i v e s ( S - S i ) ' y * I t i l ( r ) ' [ B 1 * a t 1 i ( r ) - | - B 2 * a t 1 i ( r - d ) ] d 3 r - J * i t ] 1 ( r ) A C / ( r ) [ B ^ a t . i W + B 2 * a t , i ( r - d)] d 3 r + X>lk"R ~ E^ I*lt,i(0 [5i*at,i(r - R) + B2*at,i(r - d - R)] d 3 r R#() ^ - y * i t ] 1 ( r ) A C / ( r ) [Bi* a t , i (r - R) + B 2 *at , i(r - d - R)] d 3 r J = 0 . A t t h i s p o i n t , a s e c o n d a p p r o x i m a t i o n i s i n v o k e d . N a m e l y , a s s u m e t h a t t h e w a v e f u n c t i o n d e c a y s s u f f i c i e n t l y f a s t t h a t a n y t e r m i n v o l v i n g a p r o d u c t o f w a v e f u n c t i o n s , c e n t r e d a t p o s i t i o n s t h a t d i f f e r b y m o r e t h a n t h e n e a r e s t - n e i g h b o u r s p a c i n g , v a n i s h e s . M a k i n g a c h a n g e o f v a r i a b l e s s u c h t h a t c = —d — R, a n d u t i l i z i n g t h e n e a r e s t - n e i g h b o u r a p p r o x i m a t i o n , y i e l d s {E - Ei + p} B! + j e - i k ' ( c + d ) \(E - Ex) a(c) -7 (c ) ] 1 B2 = 0 , where c 6 nn implies that the summation is over all of the nearest neighbours, ^ = - / * i t , 1 ( r ) A C / ( r ) * a t , i ( r ) d 3 r , o(c) = J * i t i l (r)* a t i l (r + c)d3r, 7(c) = J * i t i l ( r ) A [ / ( r ) * a t ! l ( r + C ) d 3 r , and it is assumed that the ^at . i are orthonormal. If the same procedure is repeated, replacing Eq. (2.2) with J ^ t ] J . ( r - d ) / f a t * ( r ) d 3 r + J **y(r-d)Ay'(r)$(r) d 3 r = J ^ . ( r - d)E9(r) d 3 r , a second equation in terms of B\ and B2 emerges. In matrix form, where the arguments from a and 7 have been dropped since their values do not depend on which neighbour we are considering, the equations may be summarized as E-Ex+P 5Ze- i M c + d>[(£--Ei)a-7]] cGnn J2 e i k-( c+d) [(E - E{)a - 7] E-E1+B .cGnn 0 B2 0 In order for this equation to have a non-trivial solution, the determinant of the 2 x 2 matrix must be set to zero. Solving yields E = Ei + V -B±1 Eeikc c€nn l±a e i k c cGnn J where the plus and minus signs go together to give the bonding and antibonding ir bands, respec-tively. For graphene, 7 = —3.033eV and a = 0.129 [4]. B is of lesser importance since it does not multiply any function of k. In addition, it is convenient to reference energies to E\ — 8, and neglect the a correction to the B term since both B and a are small. In fact, it is very common for a to be completely neglected. This is referred to as the Slater-Koster scheme [4], and this approximation will be used throughout the remainder of this work. Note that 7 and a are not, in general, exact, and considerable variation may be found in these numbers. For example, a 7 value of —2.5eV has been reported in Ref. [5]; for carbon nanotubes, Ref. [67] suggests that the appropriate value is —2.7±0.1eV in order to get good agreement between experimental bandgap measurements and the theoretical predictions. In the end, 7 is essentially a fitting parameter to the models. 2.2 Zone-Folding For a carbon nanotube, the "zone-folding" scheme is utilized to compute the dispersion, or E-k, relation from that of graphene [6]. In the Slater-Koster scheme, and referencing energies as discussed above, the tight-binding Hamiltonian may be written as H = 0 k c T E e * c £ n n o c £ n n (2.3) where the eigenvalues of this matrix, as functions of k, are the eigenenergies, E. It is convenient to work in a coordinate system where one axis is aligned with the CN axis, and the other is in the circumferential direction. If the standard rcy-plane is oriented such that the primitive translation vectors, mapping to equivalent sites in the graphene lattice, are given by Eqs. (1.1) and (1.2), then the vector L, wrapping around the circumference, is given by x 3a a\/3 . . „ L = — (ni + n 2 ) x + (m - n2) y = n ^ i + n 2a 2. Similarly, the vector T, which points to an equivalent site along the axis of the CN, including the helical angle, is given by [6] _ 3a / n 2 - ni \ „ 3a\/3 (n\ + n2 2n2 -I- n\ 2ni + n 2 y = — 3 a i J a 2 . d.R dR where dn is the greatest common divisor of 2n\ + n 2 and 2n2 + iii. It is useful to note that rotations from the xy-plane to the LT-plane may be performed via V3 27? (ni +n2) Ti! - n 2 n 2 — rij V3 27? vx VL .Vy. VT 2T? 2r? where Vi denotes the i-component of some arbitrary vector v, and 7? = y/n\ + n \ n 2 + TI^. Note that L and T span the real-space primitive unit cell of the carbon nanotube, i.e., an entire CN may be tiled, with no gaps or overlaps, by translating this cell through integer multiples of these two vectors. Since this cell is bigger than that of graphene, the tight-binding derivation detailed above would involve a larger Hamiltonian matrix than in Eq. (2.3), with one row for each atom in this enlarged cell. Here, the zone-folding scheme is employed in order to calculate the CN band structure directly from the two-dimensional (2D) result for graphene. This scheme provides the identical result to that generated from treating each atom in the CN primitive unit cell, although the full tight-binding approach lends itself to performing the simulations discussed in Chapter 5. Following the description in Ref. [4], the number of primitive graphene unit cells that fit inside the enlarged CN primitive unit cell is N = IL x T l | a j x a 2 | which means that there are 2N carbon atoms in the CN primitive unit cell. The reciprocal lattice vectors of graphene, bi and b 2 , are computed by utilizing the defining relation &i • hj = 2-ndij, where dij is the Kronecker delta function. The result is that 2 7 r ( 1 - -b i = —T= —7=X + y av /3Vv /3 i 2TT / 1 . „ Performing the analogous procedure for the reciprocal lattice vectors of the CN, KT and KL, yields K T = jj. (n2hi - nib 2), 1 / 2ni + n 2 2n2 + ni As in Ref. [4], note that N~KL is a reciprocal lattice vector of graphene, as may be seen based on the definition of (IR. For the CN, jK.i are all equivalent, where j is an integer. In addition, for j 6 {0,1,..., N — 1}, JKL are all reciprocal lattice vectors of the CN, but not of graphene. As a result, these discrete jKL in the graphene E-k relation may simply be used for the CN dispersion relation. Hence, one-dimensional (ID) projections of the graphene relation are "folded" into the first Brillouin zone of the CN. For further details, see Ref. [4]. Specifically, if k = kTT + kLt in Eq. (2.3), then where k^L = JKL and E = ±7 kT e k c c€nn 7T 7T ] x | ' i x i 10 The DOS, D(E), may then be computed numerically as Y ! TV ^ dkT 3=0 where spin degeneracy is included, and E^) is the energy associated with kiL = JKL. The E-k relation and DOS for a (16,0) CN are shown in Fig. 2.1. Each band in Fig. 2.1(a) corresponds to a different value of j, and some of these bands are multiply degenerate, i.e., different values of j produce the same E-kr relationship. In Fig. 2.1(b), note, too, the singularities at the band edges, as expected for a ID DOS. 11 9 Energy (eV) (b) Figure 2.1: (a) Dispersion relation and (b) density of states for a (16,0) carbon nanotube. 12 Chapter 3 Equi l ibr ium The simplest case to consider for C N F E T s is that of equilibrium. The goal is to compute the charge on the C N as a function of the applied potentials on the source, drain, and gate electrodes. From this relationship, the energy band diagram may be computed in order to gain some insight into C N F E T operation. A Schottky-barrier device is considered, with metallic source and drain electrodes that contact an intrinsic, semiconducting nanotube. This work is the foundation for many subsequent efforts into non-equilibrium models [51,53,55,68-76], and was initially used to study the closed [77] and open [71] cylindrical structures depicted in Fig. 3.1. The gate is translucent in order to view the inner structure of the device, and the source and drain cap either end of the inner C N cylinder. These coaxial structures are considered due to their enhanced Schottky-barrier modulation [22,78], and in order to investigate the likely performance limits of these devices, particularly with respect to their short-channel behaviour [21,23,52,70,79]. (a) (b) Figure 3.1: Schematics of the (a) closed and (b) open cylindrical C N F E T geometries. 13 3.1 Electrostatics Consider, first, the electrostatics problem. Poisson's equation is X72V = -l(p-n + Nf), where V is the local electrostatic potential, q is the magnitude of the electronic charge, p is the hole concentration, n is the electron concentration, Nf is the concentration of fixed charge, and e is the permittivity of the material. The boundary conditions, referencing applied potentials to the source, are <4>n VD = V D S - — , Q • VG = VGs-—, Q where Vs, VD and VQ are the boundary conditions on the source, drain, and gate electrodes, respec-tively, VDS a n d VQS are the applied drain-source and gate-source potential differences, respectively, and $ 5 , <!?£> and $G a r e the workfunctions of the source, drain and gate metallizations, respectively. The boundary conditions on the source and drain are appropriate in the absence of Fermi-level pinning [80]. Note that VDS = 0 in equilibrium. The fixed charge is an input to the model, and will be set to zero unless otherwise stated. The electron and hole concentrations are computed according to the model being employed, and this calculation will be detailed when each model is discussed. At this point, it is sufficient to state that the CN charge is assumed to be confined to the surface of the CN in the form, p - n = 5{p - Rc) (P2D - n 2 D ) , (3.1) where 6(-) is the Dirac delta function, p is the distance from the CN axis, and p2o and n 2 D are the two-dimensional carrier densities on the surface of the CN. 3.1.1 C l o s e d G e o m e t r y The simplest geometry to consider for the equilibrium electrostatics of a cylindrically-gated CNFET [77] is one where the source, drain and gate form a closed cylindrical box, as in Fig. 3.1(a). This case 14 is mathematically convenient since there are only Dirichlet boundary conditions bn a finite domain. In this case, it is natural to work in cylindrical coordinates with V{p,<i>,z) = V{p,<l> + 2ir,z), and finite V(0,4>, z), where <j> is the azimuthal angle, and z is the distance along the CN axis. If p 2D and ri2D are not functions of 4>, then V will also be independent of 4>. If the permittivity is constant throughout the device, it is convenient to solve Poisson's equa-tion using superposition. The first problem is that of Laplace subject to the boundary conditions described above. This solution is analytic, and is derived in Appendix A. The solution to Poisson's equation subject to homogeneous boundary conditions may be added to the Laplace result in order to obtain the full solution. The Poisson problem is solved via the Green's function, which is de-rived in Appendix B, and the electrostatics problem is reduced to one of numerically integrating the Green's function multiplied by the charge distribution on the surface of the CN. If the permittivity is not constant throughout the device, numerical methods are required. Typi-cally, it is assumed that a material of constant permittivity surrounds the CN, and that the interior of the CN itself has unity permittivity [81]. Due to the permittivity discontinuity, there exists a matching condition at the CN surface [82]: £ o x VV\R+ h-ec VV\R- • n = -qP2D " n 2 D , (3.2) c c e0 where eox and ec are the relative permittivities outside and inside the CN, respectively, eo is the permittivity of free space, RQ and RQ are infinitesimally outside and inside the CN, respectively, and n is the unit outward normal. In cylindrical coordinates, and assuming no ^ -dependence in p 2 D and n 2 o , d2V I d V d2V n is solved in two regions, p > Rc and p < Rc, subject to Eq. (3.2) across p = Rc-3.1.2 Open Geometry For the physical case, where gaps will be present between the end contacts and the gate, as shown in Fig. 3.1(b), an analytical solution cannot be derived. One option to deal with this case is to simply 15 apply homogeneous Neumann boundary conditions in the gaps sufficiently far away [83] from the surfaces of.interest, i.e., away from the CN surface. If cylindrical symmetry is present, however, an alternate solution may be employed [71]. The 2D electrostatics may be solved by applying a conformal mapping to the far-field region [84]. This yields two domains which must be solved simultaneously. The near-field region is a central disk of radius RB, which encloses the contacts and the CN, and the far-field region is the entire plane outside of this disk. Since there are no electrodes in the far-field region, the appropriate equation to solve is Eq. (3.3), where the p = RB boundary introduces a matching condition between the two domains, where it is required that V is continuous and smooth. If U = z + ip and T = £ -f ia, the conformal mapping ^_ R R T = — may be used in the far-field region in order to obtain a finite domain of the same size and shape as the near-field domain. Here, U and T represent the coordinates in the untransformed and transformed domains, respectively. Applying this transformation to Eq. (3.3) yields where the overbar indicates evaluation in the transformed coordinate system, and the derivation of this expression may be found in Appendix C. Continuity is easily satisfied on both boundaries since the transformation does not modify the value of V. The smoothness condition is given by [71] VuV • h u = -VTV • n r , where Vu and V r are the gradient operators in the untransformed and transformed domains, respectively, and and fir are the unit outward normals to the corresponding domains. By posing the open boundary electrostatics problem as two coupled closed boundary problems, Poisson's equation has been cast in a form that lends itself to the finite element technique for solution [71]. 3.2 Solution Procedure Whichever type of cylindrically-symmetric geometry is considered, y>2D and ' « 2 D must be defined. In equilibrium, this is a straightforward matter if it is assumed that V rigidly shifts the zero-bias band 16 structure as derived in Chapter 2. Due to the symmetry of the conduction and valence bands in the Slater-Koster scheme, P2D - "2D = D{E) [f {E + qV + $c) -f(E-qV- $c)] dE, (3.4) where $c is the workfunction of the CN, and f(E) = 1 exp + 1 is the Fermi-Dirac distribution function, where Ep is the Fermi level, fcg is Boltzmann's constant, TK is temperature, and E = 0 corresponds to the lowest band edge in the density of states. Initially, assume that there is no permittivity discontinuity, and consider the closed geometry. In this case, the electrostatics problem may be solved using the results derived in Appendices A and B. The solution in this case may be summarized as follows: 1. Solve Laplace's equation using Eq. (A.9) from Appendix A to get VL-2. Compute the charge via Eq. (3.4), where V = with the superscript indicating the iteration number, and = VL-3. Compute the potential due to the charge VQ\ via Eqs. (3.1) and (3.4), using a numerical integration scheme, such as Simpson's rule, to evaluate Eq. (B.2) in Appendix B. 4. Check that = VL + VQ1 is consistent with the charge to within some specified tolerance. If not, increment i, and return to Step 2. If it is, convergence has been obtained. If there is a permittivity discontinuity, or open boundaries, numerical methods are utilized. In Refs. [68,71,77], the finite element method was employed, through the use of FEMLAB [69,85], to solve the nonlinear Poisson equation that results from casting p and n in Eq. (3.1) as functions of V via Eq. (3.4). 3.3 Results and Discussion Results are presented, from Ref. [77], for the CNFET depicted schematically in Fig. 3.1(a) with a (16, 0) CN with a radius of 0.63 nm, a length of 100 nm, and a gate workfunction of 4.5 eV. The 17 electron affinity for the carbon nanotube is taken to be 4.2 eV, based on a workfunction of 4.5 eV [7], and an intrinsic-tube bandgap of around 0.6 eV. These values correspond to so-called "positive barrier" devices for electrons since the workfunctions of the source and drain are greater than the electron affinity of the CN. Similarly, "negative barrier" and "zero barrier" devices could be defined. For these positive barrier results, tc is taken to be the same as that of the gate dielectric. This permits the utilization of the analytic electrostatics solution. The temperature is taken to be 300 K. At equilibrium, if $s and $£> are equal, the energy bands along the tube will be symmetrical; thus, only profiles near one contact need to be shown. Fig. 3.2 shows these diagrams near the source for RQ/RC = 10, where RG is the gate radius, and data is presented for eox = 3.9, $s = VDS = 0, and VGS — 0-2 and 0.5 V. In Fig. 3.2(a), $ 5 = 4.5eV, which corresponds to the case of equal workfunctions for the metal and the nanotube, whereas $ 5 = 4.33eV, as in Fig. 3.2(b), and <E>s = 4.63 eV, as in Fig. 3.2(c), refer to low- and high-metal workfunctions, respectively [80]. The potential in the body of the tube, away from either of the end contacts, depends directly on VGS, leading to potential spikes near the source and drain of height determined by both 3>S,D and VGS-In the low-$s case at low VGS, thermionic emission will likely be the most significant contribution to the source current. In all other cases shown in Fig. 3.2, tunneling of electrons through the energy spike may also be important. The band diagrams for the same workfunction cases as used in Fig. 3.2, but for RQ/RC = 50, are shown in Fig. 3.3. The reduced band bending in the tube at the contacts, due to poorer coupling between the gate and the nanotube, is very evident, and will lead to a dramatic decrease in the tunneling current. Note that the present state of the art as regards gate-dielectric thinness is 2nm[50]. Regarding the permittivity of the dielectric, Ref. [86] has reported the use of zirconia, for which tox is around 5 times higher than that used in obtaining the above figures. The effect of such a change in eox can be seen by comparing Figs. 3.2 and 3.4. At VGS — 0.5 V, the increased capacitative coupling between the gate and the tube drives the mid-tube potential energy to lower values, yet does not change significantly the width of the source barrier at its base; thus, obviously, an increased current for a given bias will result from using a higher eox. At lower VGS, e9-, 0.2V, the increased eox makes essentially no difference to the band diagram because, at least for RQ/RC = 10, there is virtually no charge induced on the tube. In this case, the electrostatics are governed by the Laplace solution alone. 18 0 20 40 0 20 40 Distance from source (nm) Figure 3.2: Equilibrium energy bands near the source for a 5.6 nm gate-oxide thickness, with $ 5 = $z> set to: (a) 4.5 eV, (b) 4.33 eV, and (c) 4.63 eV. Shown are VGs = 0.2 (dashed) and 0.5 V (solid). Energies are with respect to the Fermi level (dotted) [77]. 0 20 40 Distance from source (nm) Figure 3.3: Equilibrium energy bands near the source for a 30.9 nm gate-oxide thickness, with $s = $r> set to: (a) 4.5 eV, (b) 4.33 eV, and (c) 4.63 eV. Shown are VGS = 0.2 (dashed) and 0.5 V (solid). Energies are with respect to the Fermi level (dotted) [77]. 19 0 20 40 0 20 40 Distance from source (nm) (c) 0 20 40 Figure 3.4: Equilibrium energy bands near the source for a 5.6 nm gate-oxide thickness with a per-mittivity five times higher than in Fig. 3.2, with $g = $^ set to: (a) 4.5eV, (b) 4.33eV, and (c) 4.63eV. Shown are VGS = 0.2 (dashed) and 0.5 V (solid). Energies are with respect to the Fermi level (dotted) [77]. From Figs. 3.2, 3.3 and 3.4, it appears that the width of the potential barrier at its base depends strongly on RQ for the contact geometry considered here, as has been remarked upon elsewhere [87], and here it is indicated that it also has a very weak dependence on VGS- This has been used in the development of a compact model for the non-equilibrium behaviour, where an approximate expression describes the band bending near the contacts [51,53]. Note that the effect of changing the gate workfunction from the value of 4.5 eV used here can be readily appreciated from the foregoing figures as an increase in $G of 0.1 eV, for example, has the same effect as a corresponding decrease in VGS-20 Chapter 4 The Effective-Mass Method Out of equilibrium, the modeling is complicated by the fact that Fermi distribution functions cannot be used to describe the energy spectrum of carriers in the CN. One solution is to neglect the effect of the CN charge on the electrostatics, and simply use a Laplace solution to define the non-equilibrium band structure [50]. Alternatively, one could employ infinite CN electrostatics, and consider only thermionically emitted carriers [52]. This approach is facilitated by the use of the so-called "quantum capacitance" [88] which describes the accumulation of CN charge due to changes in the local electrostatic potential [89]. Typically, this capacitance is treated as a constant [11,23, 28,30,90,91], although this is only strictly true when the charge is accumulating in a linear region of the E-k relation [89]. While this infinite-CN approach may be reasonable for certain devices, it does not allow for the treatment of carrier injection via tunneling. Tunneling has been included, via the WKB approximation, in compact models [51,53], and in a quasi-Fermi level approach for CNFETs [68,69] which follows a method similar to that used with success in HBT modeling [92]. The latter method, however, appears to predict earlier breakdown in the device saturation than the more rigorous effective-mass Schrodinger-Poisson solution presented here [72]. In this chapter, a modified WKB result is derived in order to account for the non-zero reflection of carriers above a potential energy step, and an algorithm for solving the coupled Schrodinger and Poisson equations is presented. While FEMLAB [85] had been used for some modeling of CNFETs [68,69,71,77], the program is somewhat unwieldy when making geometrical changes to the device structure, or when performing batch simulations. For these reasons, this chapter begins with a description of a finite difference solution for the Poisson equation. 21 4.1 Finite Difference Electrostatics The finite difference method for solving Poisson's equation is well-established in Cartesian coordi-nates. For a closed cylindrical geometry, however, some modifications are required. Recall that Eq. (3.3) is being solved in two regions, subject to the matching condition given by Eq. (3.2) across the CN surface. The finite difference approximation to Eq. (3.3) is Vj+i,j ~ 2Vj,j + Vj-i,j , Vj+i,j - Vi-i,j , Vj,j+i - 2Vj,j + Vjj-i _ n i \ Ap2 + 2piAp + Az2 _ U ' 1 1 where Ap and Az are the uniform mesh spacings in the p- and 2-directions, respectively, and the subscripts i, j give the node numbers in the radial and longitudinal directions, respectively. Note that, since the mesh spacing is assumed uniform, Pi may be written as iAp, and that 0 < p < RG-This equation can be solved for Vij for nodes where p ^ 0 and where p ^ Rc-For iAp = Rc, the finite difference approximation to Eq. (3.2) is eox {Yi±l^_Yi^ _ £c (ViJ = - g P 2 D % " Q n 2 D J , (4.2) where p2D,j and ri2u,j are the 2D hole and electron concentrations at the jth node in the z-direction on the CN surface. For i = 0, i.e., on the CN axis, symmetry implies that 8V dp = 0. p=0 In this case, l'Hopital's rule is applied to the singular term of Eq. (3.3). This yields l i m ( - — ) = — p—*o \p dp ) dp2 p=0 For i = 0, then, the appropriate equation is Ap2 + Az2 ' 1 1 where V-ij = Vij by symmetry. Eqs. (4.1), (4.2), and (4.3), may be solved simultaneously by formulating the system of equations in matrix form. Note that it may be helpful to normalize these equations such that the coefficients of the matrix do not vary too greatly in magnitude. This is to aid the numerical matrix inversion. 22 4.2 Treating Quantum Phenomena With the Poisson equation solved, consider, now, the computation of the CN charge. In equilibrium, this was facilitated by assuming that Fermi-Dirac statistics applied throughout the tube. Out of equilibrium, this is not the case. If the charge were either known or negligible, one could solve the electrostatics problem, and then simply compute the current corresponding to that band diagram. If ballistic transport is assumed, the current X is given by the Landauer equation [93,94] where T(E) is the transmission probability for the scattering region, i.e., the channel, h is Planck's then, it is sufficient to compute T{E) in order to determine the terminal current. 4.2.1 A s y m p t o t i c Treatment of Tunne l ing For the tunneling problem, the WKB approximation is perhaps the most widely used approach. The method is derived in Appendix D, with the standard result that where U{z) is the band edge, h is Dirac's constant, mcff is the carrier effective mass, and the turning points, z0 and z\, are such that U{z0) = U{z\) = E with U(z) > E for z0 < z < z\. If there are no points where this occurs, and E > U(z), then T(E) = 1. If there is only one point where this occurs, or E < U(z) everywhere, then T{E) = 0. Note that this expression assumes that U(z) is slowly varying. If, as a further approximation, Eq. (4.5) is taken to be an equality, a simple T(E) is obtained that may be used to compute the terminal current [51,68,77]. While this approximation may be useful, it can be shown, at least for Schottky barrier devices where there may be an abrupt change in U(z) at the contacts, that this may result in a significant overestimate of the current [53]. This is due to the possibility of reflection from the potential energy step even for carriers above the barrier. These carriers would be assigned a unity T(E) in the above scheme. A better expression for this thermionic current component may be derived while still using the WKB ansatz [53]. (4.4) constant, and spin- and CN band-degeneracy have been included. If the band diagram is known, (4.5) 23 For phase-incoherent transport, the transmission through the regions close to the source and drain contacts may be computed separately. The total transmission probability may then be computed from T*{E) TSTD TS + TD- TSTD' where Ts and T/j are the source and drain transmission probabilities, respectively, and the star superscript indicates that phase-incoherent transport is being considered. For specificity, the WKB transmission probability is derived for a source electron with energy E greater than any point in a "barrier" region of width w, where U(z) has a significant discontinuity at some point. In particular, consider the case of the source barrier, as in Ref. [70], which has a band • diagram that is similar in form to that shown in Fig. 4.1. Note that this figure is intended to be generic, and that the exact same calculation will hold for a barrier that has the opposite concavity, such as in a positive-barrier device, or for a barrier profile that is less smooth, as long as it is smooth enough for WKB to hold in the barrier region itself. w Figure 4.1: Model energy band diagram. There are several important features of this energy band diagram: a discontinuity occurs at z = 0, the potential energy is constant for z < 0 and z > w, the potential energy and its derivative are continuous across z = w, carriers are only incident from the left, and one effective mass describes the carriers in all three regions. 24 The usual WKB form for the solution can be written as [95] Ci exp (iksz) + C2 exp (-iksz) , 1 , x C 3 e x p ( i / kB(z) dz ) + C4exp ( - i / fcB(z)di \fkB(z) V \ Jo J V Jo C5exp [ i fcM (2 - w)] z < 0, 0 < 2 < W , 2 > I f , where fcs, fcs, and fcjw are the wavevectors when z < 0, 0 < z < w, and z > w, respectively, and Ci through C 5 are constants. Note that kB(w) = kM and k'B(w) = 0 from the key assumptions of this model problem, where the prime denotes a total derivative. In this case, only the transmission probabilities are of interest, so C 5 may be arbitrarily assigned a value of unify. The usual matching condition [96] is that ^ and its derivative must be continuous. The deriva-tives of the wavefunction are iks [Ci exp (iksz) - C2 exp (-iksz)] C3e" 2ik2B 2k j C4e-2ikB -\- kB 3 2ka ikMC5 exp [\kM(z - w)] z < 0, 0 < z < w , z > w, where Idi, v= j kB(z)i Jo and the explicit z-dependencies of kB, k'B, and v have been dropped for notational convenience. Note that z is a dummy integration variable used to avoid confusion with z in the integral limit. These expressions may be further simplified. Since the barrier region joins smoothly onto the right-most region, it is expected that there will, in fact, be no reflected wave anywhere in the barrier region; the initial discontinuity at z = 0 is the only source of reflection. This essentially stems from the form of the solution assumed for the output wave, and from the fundamental WKB assumption. Of course, the same result will be arrived at even if C4 is not set to zero at this point; there will just be a more difficult set of equations to solve. The simplified system of matching equations is C3 Ci + c2 = C 3e i l / ( l u ' fk^ = C 5 , 25 i b ( C l - C ! ) = C l M V 24(0) -with the solution, upon setting C 5 = 1, Ci = [2kskB(0) + 24(0) + ifc^O)] , 4kskB(0) C 2 = [2kskB(0) - 2k2B(0) - ifcj,(0)] , 4kskl(0) C3 = V ^ e - ^ ) . The source transmission probability is given by ks In this case, therefore C l ' 2 = 16fcl5(0) { ^ ( O ) ] 2 + 4 [fc|(0) + feMO)]2} , 16fcs4(0) [^(0)]2+4[fc|(0) + fcsfcB(0)]2 Ts = , ,T* . . , „2. (4-6) which reduces to unity if k'B(0) is sufficiently small, and kB(0) = ks, i-e., when U(z) has no discontinuity, and is slowly varying. An analogous expression holds for TQ-Incorporation of reflections from the abrupt potential energy step, while still using WKB, results in much better agreement with a full Schrodinger solution [53], and is much less computationally expensive than the full solution. 4.2.2 Solving Schrodinger's Equation In order to perform a full solution to Schrodinger's equation [72], the scattering matrix method [97] is employed, where input wavefunctions are assumed at either end of the scattering region. The transmission matrix method [98] has also been attempted, where an input wavefunction is assumed at one end, and an output wavefunction at the other. This latter approach, however, while conceptually and implementationally simpler than the former, becomes numerically unstable for low transmission probability wavefunctions. As a result, the scattering matrix method, which has been sufficiently described in Ref. [97], is preferable. 26 In this system, the charge density is given by O I D -n1D)S(p-Rc) • P ~ n = 2^p- ' where nm(z) and P I D ( Z ) are computed via <92* 2m e f f dz2 h? (E - 17) ¥, (4.7) where the nanotube effective mass is obtained by fitting a parabola to the tight-binding approxi-mation of the band structure, and is the same for both electrons and holes due to symmetry [4]. The fit for the lowest conduction band of a (16,0) C N is depicted by the dashed line in Fig. 4.2. This approximation is reasonable for energies near the band minimum. For higher energies, the Figure 4.2: A parabolic fit (dashed) to the (16,0) C N band structure (solid). non-parabolicity of the bands is clear from the figure. Only the first, doubly-degenerate band is included in the calculations presented herein. The potential energy, or band edge U, for each carrier type is specified by Ue(z) =Emc{z)-xc, Uh(z) = -Ue(z)+EG, 27 where EG and \c are, respectively, the nanotube bandgap and electron affinity, EVSbC is the vacuum level which is given by — qV, and the e and h subscripts refer to electrons and holes, respectively. Note that the sign on Uh is chosen such that a subsequent integration will be over positive energies in order to get the hole concentration in the valence band. The scattering matrix method provides a numerical solution for Eq. (4.7) by cascading 2x2 matrices [97]. In the required discretization of the potential energy, it is found that the use of a piecewise constant approximation, with plane-wave solutions, is preferable to a piecewise linear approximation, with Airy function solutions, due to the considerable reduction of simulation time without an appreciable increase in the error. Matching of the wavefunction and its derivative on the boundary between intervals i and i + l, assuming a constant effective mass, is performed via the usual continuity and smoothness conditions. In order to completely specify the wavefunction, two boundary conditions are required. In the contacts, the wavefunction at a given energy is of the form where ks and ko are the wavevectors in the source and drain contacts, respectively, Lc is the CN length, as depicted in Fig. 3.1(a), and C\ through C4 are constants. As an example, noting that an analogous calculation may be performed for the drain by exchange of variables, the case of source injection is now illustrated. For this case, C 4 = 0 for all energies. This is one boundary condition. In addition, it is expected that the Landauer equation [93] will hold for the flux, and must be equal to the probability current. For the transmitted wave, this yields [72] where the pre-factor of 2 accounts for the aforementioned band degeneracy, fs is the Fermi-Dirac carrier distribution in the source, and T is specified by C i e i f c s 2 + C2e-[ C3eikDZ + de'' *LfsT=*LkD\ch\\ T = kp\C3\2 Simple manipulation yields the normalization condition 2meff fs ith2 ks (4.8) 28 At each energy, multiplication of the unnormalized wavefunction by a constant satisfies Eq. (4.8). Including source and drain injection components, the normalized wavefunctions yield the total carrier densities in the system, n1D(z)= ^ (|* e,s| 2 + | * e , D | 2 ) dE, Pm(z) = / (l*fc,s|2 + |*fc,D|2) dE, J £h where the subscripts e and h refer to conduction and valence bands, respectively, S and D refer to source and drain injection, respectively, and £e>h is taken to be the bottom of the band in the appropriate metallic contact. In this way, penetration of the carrier wavefunctions from the metal into the bandgap of the CN can be included [99]. In practice, the integrals are performed using adaptive Romberg integration, where repeated Richardson extrapolations are performed until a predefined tolerance is reached [100]. An adaptive integration method is a necessity in order to properly capture ty, which is typically highly peaked in energy for propagating modes. Without adequate sampling of the wavefunction, convergence is unlikely. Alternatively, one could employ a very fine discretization in energy; however, the Romberg method allows for the mesh size to change based on the requirements of the integrand, and results in a much improved simulation time. As a final note to this section, a check is performed to confirm that the normalization results in the expected form for the electron and hole concentrations when there is no band bending anywhere. In this case, fs = fo = f(E), where fr> is the Fermi-Dirac distribution function in the drain; ks = ko = k(E), and T = 1. Considering electrons, n i D { z ) = k ^W) Substituting in the definition for k(E), i.e., yields rOO n1D(z)= / 2D(E)f(E)dE, j£e where D(E) is the standard ID, single band, effective-mass DOS, and the factor of two accounts for the double degeneracy of the lowest conduction band. This is indeed the expected expression. 29 4.2.3 Quantum Capacitance Since a self-consistent charge-voltage relationship is being computed, it is instructive to consider the "capacitance" of a carbon nanotube [89], i.e., the sensitivity of the charge accumulation to the local electrostatic potential. The term "quantum capacitance" has been used in the CNFET literature [30,90] in order to encapsulate this behaviour, and was most notably used by Luryi [88] in order to develop an equivalent circuit model for devices that incorporate a highly conducting two-dimensional electron gas. In the CNFET case, the results are somewhat different [89]. Here, this potentially important concept is discussed in the context of CNFETs, including a derivation of the 2D case, which has been previously presented in Ref. [88], in order to help illustrate some key differences. In order to derive analytical expressions, it is assumed that the CNFET is in quasi-equilibrium, and that the carrier distribution functions, as usual, are rigidly shifted by the local electrostatic potential. The quasi-equilibrium assumption permits the use of Fermi-Dirac statistics for the carrier distributions in energy, and is reasonable for the case when tunneling is not significant. Since the DOS is symmetric with respect to EF, the net carrier density, due to electrons and holes in the semiconductor, may be written as rOC P d ~ n d = / D{E) Jo where the subscript d represents the dimensionality of a quantity, Vcs is the local electrostatic potential, EG is the bandgap, and Ep is taken to be mid-gap when Vcs = 0. The quantum capacitance, CQ, is defined as ° Q = q dVcs ' ( 4 - 1 0 ) and has units of F/m 2 and F/m in the 2D and ID cases, respectively. In the 2D case, if the effective-mass approximation with parabolic bands is employed, the DOS is given by °{E) = SM(S)' where M(E) is the number of contributing bands at a given energy, and the meff for the conduction and valence bands are assumed to be equal as in ID CNFETs. If this is combined with Eqs. (4.9) E + ^ + q V c s ) - / E + ^ - q V c s dE, (4.9) 30 and (4.10), and the order of differentiation and integration is exchanged, the result is that CQ = metrq' - / M(E) • K Jo sech 4Trh2kBTR If M is a constant, the integration yields E + ^ - q V C S \ + g e c h 2 (E+fy+qVcs 2kBTA Cc = Mmenq2 2nh2 K sinh 2kBTK dE. 2 -\2kATK J cosh 1kBTK 1 CUB" 2kBTK J which reduces to CQ = MmeSq2 •Kh2 (4.11) when EQ = 0, in agreement with Ref. [88], where metallic properties were assumed. Note that this function is piecewise constant in the metallic case, i.e., it is quantized. In the semiconducting case, the function has continuous variation in it for a given M. For EQ greater than about 15kBTx, however, the function makes a rapid transition from a small value to that given by Eq. (4.11) when Vcs crosses EQ/2, and is thus effectively quantized. In the ID, effective-mass case, D(E) = M{E) /2mef f The explicit energy dependence of this DOS complicates the evaluation of the integral for CQ. The approach suggested in Ref. [30], using the fact that the derivative of f{E) is peaked about EF in order to approximate this integral using a Sommerfeld expansion [66], cannot be followed in general, due to the presence of singularities in the ID DOS. This approximation may be useful if attention is restricted to situations where the DOS singularities are sufficiently far away from EF-The capacitance is given by CQ Iruefi 2kBTKhM 2 f Jo M(E) s/E sech +sech2iE+ 2 (4.12) For sufficiently large \Vcs\, one of the sech2(-) terms may be completely neglected. As a simple example, if Vcs = 0.1 V for a material with Ec — leV, the contribution to the integral from 31 the first term is roughly four orders of magnitude greater than the second. This approximation is equivalent to neglecting hole charge for positive Vcs, and electron charge for negative Vcs- The solid line in Fig. 4.3 shows the equilibrium CQ as a function of Vcs for an arbitrary ID semiconductor with two valence and conduction bands: at 0.2 and 0.6 eV away from the Fermi level. An effective mass of 0.06mo is assumed, where mo is the free-electron mass. The van Hove singularities, at each band edge, result in corresponding peaks in CQ. Figure 4.3: ID quantum capacitance as a function of the local electrostatic potential at equilibrium (solid), and for the mid-length region of an end-contacted semiconductor with a bias voltage of 0.2 V (dashed) between the end contacts. The effective mass is taken to be 0.06mo, and energy bands are situated at 0.2 and 0.6 eV on either side of the Fermi level [89]. For a linear energy-wavevector relation, such as that near the Fermi level in graphene or a metallic CN, the DOS is constant. This is the case considered in Ref. [90] and is valid when Vcs is such that f(E) is approximately zero before the first van Hove singularity is encountered in the integral. Since the higher energy bands are not relevant to the integration under such a condition, M is constant, 32 and the DOS is given by where vp is the Fermi velocity. The result is which agrees with the expression quoted in Ref. [90]. Note that in. Eq. (4.12), CQ does not manifest itself as a multiple of some discrete amount, so "quantum capacitance" is not an appropriate description, unlike in the metallic 2D and ID metallic-CN cases, where the capacitance is truly quantized. The discussion can now be extended to include the non-equilibrium behaviour for a general, ID, intrinsic semiconductor. All of the numerical results are based on the methods described in Refs. [70] and [72], which consider the cases when transport in the ID semiconductor is either incoherent or coherent, respectively. While these methods were developed in order to describe CNFETs, their use of the effective-mass approximation allows them to be used for any device, and bias, where the semiconductor is described well by this approximation. For phase-incoherent transport, a flux-balancing approach [51,70] is used to describe the charge in a CNFET. Considering only the electrons that are far away from the contacts, i.e., in the niid-length region, Eq. (4.9) becomes [51] where Vcs is evaluated in the mid-length region. A similar expression holds for holes. The first term in Eq. (4.13) resembles the equilibrium case, so a similar form for that contribution to CQ is expected. The peak for each contributing band will occur at the same Vcs, but the overall magnitude will be smaller due to the multiplication by the transmission function. The second term is also similar except that these peaks will now be shifted by VDS- This is depicted by the dashed curve in Fig. 4.3, where the case illustrated by the solid curve has been driven from equilibrium by VDS = 0.2 V. Note the splitting of each large peak into two smaller peaks: one at the same point, and the other shifted by VDS- Of course, the numerical value of the non-equilibrium capacitance depends on the exact geometry considered, as it will influence both Vcs a r i d the transmission probabilities in Eq. (4.13), but the trends shown here are general and geometry independent. Pd-nd= - -(4.13) 33 For the coherent, non-equilibrium case, it is instructive to consider a metal-contacted device, in which the band discontinuities at the metal-semiconductor interfaces are sufficient to allow significant quantum-mechanical reflection of carriers even above the barrier. Further, this discussion is restricted to short devices since the importance of coherence effects is diminished as the device length is increased. Due to the phase coherence, then, the structure is very much like a quantum well, even for devices where tunneling through the contact barriers is not important. For this device, quasi-bound states are expected to emerge at the approximate energies ^ i2^2h2 1 ~ 2mef[Lc' For m0ff ~ 0.06mo, such as in a (16,0) CN, Ei ~ 6.3(i/L'c)2 eV, where L'c is in nanometres. This may be compared with the result for metallic CNs, where the linear energy-wavevector relationship yields a l/Lc dependence [90,101]. Fig. 4.4(a) displays CQ as a function of position and Vcs f° r this choice of meff. The maxima, indicated as brighter patches, show a dependence on Vcs that reveals the population of quasi-bound states. Moreover, the maxima in position clearly show the characteristic modes expected from the simple square-well analogy. Note that the peak splitting occurs for coherent transport as well, as shown in Fig. 4.4(b), where the peaks have been split by Vbs= 0.1 V. The main difference, between the coherent and incoherent cases, is the presence of the quasi-bound states. These serve to increase the number of CQ peaks, since each quasi-bound state behaves like an energy band, and they also give rise to a strong spatial dependence. While Fig. 4.4 shows only a single-band, coherent result, inclusion of multiple bands would cause CQ to exhibit peaks corresponding to each band, and to each quasi-bound state. 4.3 Planar Devices Up to this point, attention has been paid only to cylindrically-gated devices. Experimentally, of course, planar devices, such as the one depicted in Fig. 4.5, are much more feasible. Without cylindrical symmetry, a three-dimensional Poisson equation must be solved for the local electrostatic potential. 34 - 0 . 1 5 -0 .15 (a) (b) -0 .2 • -0.2 -0.251 I -0 .25 > 0 10 2 0 0 10 2 0 D i s t a n c e f rom s o u r c e (nm) D i s t a n c e f rom s o u r c e (nm) Figure 4.4: ID quantum capacitance, in arbitrary units, for a short-channel, phase-coherent semi-conductor as a function of position and the local electrostatic potential for applied bias voltages of (a) 0 and (b) 0.1 V between the end contacts. The bright areas indicate higher capacitance. 4.3.1 Finite Element Electrostatics In order to simulate non-coaxial device structures, a numerical method for solving Poisson's equation, that can handle a variety of domain shapes, is desired. For a rectangular domain, the finite difference method described previously is extremely simple to implement. However, if the domain has a slightly more complicated structure, it quickly becomes unwieldy. In this case, the finite element method is chosen. Since the charge is confined to the surface of the C N , by assumption, it is sufficient to solve Laplace's equation in two regions subject to a matching condition, as described previously. If all lengths are normalized to Af, the problem to solve is •^ 2 V • (eoxVV) = 0 outside the C N , - ^ V • (ecVV) = 0 inside the C N , 35 Figure 4.5: Planar geometry CNFET with cylindrical source and drain electrodes, and planar gate electrode. The volume outside of the CN and electrodes is filled with a homogeneous dielectric. subject to the matching condition given by Eq. (3.2), except that V is divided by N to account for the length normalization. This can be written in the "weak form" by multiplying through by a test function v, and integrating over the volume of each region. The restriction on v is that it must vanish on any boundaries where a Dirichlet boundary condition applies. The equations to solve become where fl+ and il~ are the volumes outside and inside the CN, respectively. Using a standard vector identity and the divergence theorem yields -JP ( U L w w • v"d3- - I L w v " • *) - °' where dQ,c is the CN surface, d2r is a differential area element, and the other surface integrals have vanished due to the restriction on the function v. Adding these together, and rearranging, gives /fb+ e ° x W • V v d 3 r + H L e c V V ' V w d 3 r = ~ ILnv ( e ° x w | * s " e c w | * c ) • A d V -Applying Eq. (3.2) yields / / / e o x W • Vt)d3r + / / / ecW -Vvd3r = q/V [f v(P™ ~ n™) d 2 r . (4.14) JJJn+ JJJii- JJanc V £ o / 36 This is known as the "weak form" of the electrostatics problem. In the finite element method, the weak formulation is solved, and in the Galerkin technique, v and V are expanded in terms of basis functions, i.e., v = Y, V'^> i i where B{ and Vi are constant coefficients to the basis function W j . These expansions can be sub-stituted into Eq. (4.14), and since the expansion must be true for any choice of v that vanishes on the Dirichlet boundaries, the sum over the Bi coefficients may be treated term-by-term. This allows those coefficients to cancel out of the final expression, which is Vj {^Jjj^ e o xVW j- • Vwi d3r + J jecVujj • Vw; d3r j = ^ jj^ wt (P2D - n2D) d2r. This may be written compactly as ^ ] Kij Xj — Cj, . 3 or, in matrix form, KV = C, where K is known as the stiffness matrix. It is convenient that K only needs to be computed once per geometry. The forcing term, C, needs to be computed at each iteration. In order to perform these calculations, the choice of basis function must be specified. In this work, piecewise linear basis functions are selected, where Wj is unity at the i t h node, decreases linearly to vanish at the neighbouring nodes, and remains zero throughout the remainder of the domain. This is most easily seen in ID, as shown in Fig. 4.6. Note that Vi will equal the potential at the i t h node. The basis function for a 3D domain is a straightforward generalization. In this work, the freeware program TetGen [102] is used to mesh the 3D domains. 4.4 S o l u t i o n P r o c e d u r e The solution procedure is similar to that described in Chapter 3 except that the charge is now computed via a Schrodinger solution. Whether the finite difference or the finite element method is 37 i-2 1-1 I 1+1 1+2 Node number Figure 4.6: ID piecewise linear basis function for use in the finite element method. used to achieve the Poisson solution, the simulation procedure begins by precalculating the matrix corresponding to the discrete equivalent of the differential operation, e.g., the stiffness matrix in the finite element case. The solution procedure proceeds as follows: 1. Solve Laplace's equation to get Vj,. 2. Solve Schrodinger's equation over a sufficiently large energy range, and integrate to obtain the charge density. Here, EveiC = —qV^\ where V ' 0 ' = VL- In practice, this method treats one conduction band considering electrons, and one valence band considering holes. For electrons, the energy range is from the bottom of the band in the metal lead to around lbksTx above the Fermi level of the injecting contact. Source and drain injection are treated independently. As mentioned previously, the integration over energy is done using adaptive Romberg integration [100]. 3. If using the finite element method, the charge density will need to be multiplied by the basis functions and integrated in order to obtain the C vector. In the finite difference case, this integration is not required, although the charge vector still requires generation when solving the system using matrices. 4. Compute the residual r, which, in the finite element notation, is given by = KV^ — C^. 38 5. Update the potential according to y( J + 1 ) = V"W - •dK~1r^i\ where •d is a damping factor that is required in order to get convergence. A value of 0.4 was typically used for this, parameter. In practice, K^1 is not explicitly calculated, but rather, the backslash operator in Matlab is used in order to compute K~xr^ using Gaussian elimination. 6. Check that | |if" _ 1r^ ,)| | 0 0 ' is below some tolerance, where || • | |oo is the infinity norm, i.e., the maximum value of \K~1r^\. If not, increment i, and return to Step 2. If it is, convergence has been obtained, and the current may be computed via Eq. (4.4), where T(E) is stored from the computation in Step 2. Note that this iteration scheme is similar to that used in Chapter 3 when -d = 1. When t? ^ 1, this is referred to as Picard iteration. 4.5 Results and Discussion This method has been used to predict the performance of cylindrically-symmetric devices, and also to suggest the possible performance of other, asymmetric, devices. Here, results are presented for the DC behaviour of CNFETs, and a metric for the switching time, TJ, is examined, is an important figure-of-merit for digital applications, and may be evaluated using DC quantities. Small-signal results have also been generated from these DC simulations [73-76,103], and the analyses of these can be found in the appropriate references. 4.5.1 DC For the DC characteristics, as in Chapter 3, the closed, cylindrically-symmetric device shown in Fig. 3.1(a) continues to serve as the model geometry, but now ec is. allowed to differ from e o x . In particular, ec = 1 [81] for all devices, and e Q x is chosen as appropriate for the choice of dielectric. Due to the permittivity discontinuity, the finite difference Poisson solution is utilized. The model device includes a 20nm-long, (16,0) CN, and a 2.5nm-thick oxide. Additionally, all of the workfunctions are 4.5 eV, e o x = 25, and £ e , / i is taken to be 5.5 eV away from the contact Fermi level unless otherwise specified. This is equivalent to choosing a Fermi wave vector in the metallic contacts of around 0.3 A - 1 when the effective mass is equal to that of the (16,0) CN. This choice is similar to the smallest value used in Ref. [54]. 39 Beginning with band diagrams, Fig. 4.7 shows the results for Vos = 0 and 0.4 V with VGS = 0.5 V. In equilibrium, the carrier concentrations, particularly away from the source and drain contacts, are •0.2\ 0 5 10 15 2 0 D i s t a n c e f rom s o u r c e (nm) Figure 4.7: Conduction band edges for VGs = 0.5 V, and VDS = 0 (dashed) and 0.4 V (solid) [72]. in reasonable agreement with the results from using Fermi-Dirac statistics. It is interesting to note that the agreement is not exact, and this is attributed to the rapid potential energy variation at the end contacts, which leads to resonance in the carrier wavefunctions. Strictly speaking, then, a Schrodinger-Poisson solution is required even in equilibrium. Application of a drain voltage serves to decrease the number of electrons in the channel, and hence, the conduction band in the mid-tube region drops to a lower value. In both cases, due to the short CN length and large band discontinuities at the contacts, in-terference effects and contact reflections modify the carrier distribution functions in the CN. Note that this is significant even for carriers above the contact barriers. To see this, consider Fig. 4.8, where the conduction band edges and transmission probabilities for electrons with VGS = 0.5 V and VDS — 0.4 V are shown. In (a) is a "negative barrier" device, where all of the carriers in the mid-tube region are thermionically emitted. Clearly, approximating the transmission probability as unity for these carriers would be grossly inaccurate. Use of Eq. (4.6) results in a much more accurate value despite the fact that it neglects interference effects, and only includes contact reflections [53]. The same result, but for a "positive barrier" device, is shown in (b). Fig. 4.9 shows the effect of 40 0 10 2 0 0 0 .5 1 z (nm) P robab i l i t y Figure 4.8: Conduction band edges and transmission probabilities for VGS = 0.5 V, VDS — 0.4 V, with (a) $ S = $ D = 3.9 eV, and (b) $ s = $ D = 4.5 eV [72]. the resonant peaks in transmission on the carrier distribution in energy for the same device as in Fig. 4.8(a). Note that this function bears little similarity to a Fermi-Dirac distribution function multiplied by the density of states. As a result, quasi-Fermi level techniques [68,69], as used in HBT modeling [92], are difficult to rigorously justify. Consider, now, the terminal characteristics. Fig. 4.10 shows the TD-VGS characteristic for VDS = 0.4 V with <&s = $£) set to 3.9, 4.3, and 4.5 eV. Note that <3>G is left fixed at 4.5 eV since changes in this parameter, as noted previously, are equivalent to changes in VGS- An important thing tb note is that these CNFETs are ambipolar [104-107], i.e., the current is either electron- or hole-dominated depending on VQS- This is depicted schematically in Fig. 4.11, where low VGS yields a hole-dominated current, and high VGS yields an electron-dominated current. For the case shown in Fig. 4.11, where all workfunctions are equal, it is expected that the current minimum will occur at VGS — VDS/2- This is confirmed by Fig. 4.10. When the contact workfunctions are adjusted, this minimum moves, and may also change somewhat in magnitude. 41 0 . 2 0.4 0 . 6 0 . 8 N o r m a l i z e d n ( E ) Figure 4.9: Carrier concentration as a function of energy, normalized to its maximum value, for VGS = 0.5 V, VDS = 0.4 V, with $ s = $ D = 3.9 eV [72]. Figure 4.10: XD-VQS for a closed cylindrical device with Vos = 0.4 V, and = ®D set to 3.9 (solid), 4.3 (dashed), and 4.5eV (dotted). 42 0.4 0.2 > 0 -0.2 CD C L U -0.4 -0.6 -0.8 10 z (nm) 15 20 Figure 4.11: Schematic showing the ambipolar conduction mechanism. Shown are the hole (dotted) and electron (dashed) dominated cases along with the crossover (solid) between these two cases [70]. 4.5.2 Switching Speed In order to examine the suitability of CNFETs for digital applications, consider TJ, as analyzed in Ref. [55]. A method for computing this quantity from DC data, for devices with non-optimized threshold voltage [108], has been proposed in Ref. [30]. Here, this metric is examined, and results generated from the effective-mass Schrodinger-Poisson solver are presented. In order to compute r ,^ consider the circuit shown in Fig. 4.12, where a load (transistor B) is discharged by a driver (transistor A) at a rate given by the current I, defined as d Q g ( V & , V & ) di where VQS V§S = —VDD> where VDD is the positive supply voltage given by the difference between the OFF and ON gate voltages on the load, namely: VDD = VGS,OFF ~ VGS,ON-The discharge current X = —ID w ^ v a r v ^ t n e load discharges. As an approximation, the practice is to evaluate I at its highest value, namely: Z = -ID(~VDD, VGS,ON) = O N -43 DD B DD GS.OFF GS.ON Figure 4.12: Circuit for computing the delay time [55]. Expressing the change in nanotube charge in terms of a capacitance, ^ = -c8(vgs,-vDD), where Qc is the charge on the C N channel, which is equal and opposite to that on the gate if inter-electrode capacitances are neglected, i.e., for the case of the "intrinsic" switching time. As Fig. 4.13 shows, CG is not constant during the discharge period. Again, an approximation can be made by evaluating it at its highest value, shown by the slope of the dashed line, namely: CG — CG(VGS,ON> ~VDD) = CQ, ON-Now, overestimates have been made on both CG and T. Unfortunately, it is not yet clear how these approximations might compensate each other in the resulting expression for the switching speed: Td = CG,O^DD 1 A ON I The goal, of course, is to have a short switching time, and a high ratio of ON-current to OFF-current. In order to simultaneously view these parameters, the convention has been to plot TJ, against 44 Figure 4.13: Computation of the gate capacitance in the ON-state via the slope of the dashed line [55]. ION / Z O F F , where the latter quantity is evaluated by sliding a window of width VDD along a plot of XD VS. VGS, as described previously in Ref. [109]. For these results, the planar geometry depicted in Fig. 4.5 is chosen, and the finite element method is used to solve Poisson's equation. Since the overall solution method uses a ID Schrodinger equation, the axial potential of the CN is used when computing EVAC. Three devices are considered: a (22,0) CN with an 8nm oxide, a (10,0) CN with an 8nm oxide, and a (10,0) CN with a 2nm oxide. Parameters that are common between all three devices are: eox = 16, as appropriate for hafnia, the end contacts are Pd, the underlap of the gate and end contacts, as depicted in Fig. 4.14, is LSG = LDG = 4nm, the gate length, thickness and breadth is 20, 5, and 10 nm, respectively, and the end contacts have a radius i?Cont = 6nm and length L c o n t = 1 nm. The (22,0), 8nm device was chosen due to its similarity to a high-DC-performance experimental device from Ref. [110], and the other two devices are to show the effect of varying the tube diameter, and hence the bandgap, and the oxide thickness. In all cases, the workfunction of the nanotube is 4.7eV [111], and the workfunction of the gate is adjusted to give reasonable I O N in the simulation range of —1 < VQS < —0.5 V. The barrier height for the Pd/nanotube end contacts is taken from Ref. [112], and has values of —0.04 and 45 Figure 4.14: Cross-section of the planar geometry CNFET [55]. 0.3 eV for holes in the (22,0) and (10,0) cases, respectively. Unless otherwise stated, simulations are performed at VDS — —0.5 V, and 2ON/2OFF is evaluated at VDD = 0.5 V. The choice of this value is arbitrary [108], and is used here as it may be realistic for far-future, low-voltage logic applications. Fig. 4.15 shows the XD-VQS characteristics for the three model devices. While the (22,0) device has the best ON-current, due to its smaller bandgap, it also has the worst OFF-current. This is due to ambipolar conduction [105] setting in before a low OFF-current can be attained [108]. In order to suppress this effect, a larger bandgap CN is required [109] such as in the (10,0), 8nm device, where there is roughly twice the bandgap of the (22,0) CN. Due to the higher barrier height for this device [112], a significant reduction in the ON-current occurs. This may be compensated for by thinning the oxide down to 2 nm, where reasonable ON and OFF currents are observed. Fig. 4.16 displays the intrinsic switching time against the ON/OFF ratio. Due to the high ON-current in the (22,0) device, there is an extremely short switching time. Unfortunately, this device has the worst ON/OFF ratio of the three devices. By choosing the smaller CN, the ON-current is sacrificed, but a gain in ON/OFF ratio is obtained. Note that these values are in rough agreement with the sub-picosecond switching times predicted for CNFETs out to XON/^OFF * 104 in Ref. [109] and « 102 in Ref. [108]. In all three cases, although it is not apparent in the (10,0), thick-insulator case until longer values of than are shown in Fig. 4.16, T O N A O F F starts to decrease after reaching its maximum value. This corresponds to VGS.OFF passing through, in a positive direction for the p-type devices considered here, the point at which | Z D | is a minimum [108]. Additionally, the (10,0), thin-insulator case shows an increase in TJ, at the low end of its 2QN/2OFF range, which corresponds 46 " I . . . I I -1 -0.8 -0.6 -0.4 -0.2 0 F i g u r e 4 . 1 5 : XD'-VGS f o r t h r e e p l a n a r d e v i c e s [55]. t o o p e r a t i o n a t h i g h , n e g a t i v e v a l u e s o f VGS- TO e x a m i n e t h i s p h e n o m e n o n , c o n s i d e r F i g . 4 . 1 7 , w h i c h s h o w s t h e VGS d e p e n d e n c e o f t h e c h a r g e i n t h e n a n o t u b e , f r o m w h i c h t h e i n t r i n s i c c a p a c i t a n c e i s c o m p u t e d . I t i s c l e a r t h a t t h e r e i s a p r o n o u n c e d i n c r e a s e i n CG a s t h e n e g a t i v e g a t e b i a s r e a c h e s h i g h v a l u e s , a n d t h i s i s r e s p o n s i b l e f o r t h e i n c r e a s e i n T^. T h e r i s e i n c h a r g e c o m e s a b o u t b e c a u s e t h e t h i n n i n g o f t h e p o t e n t i a l b a r r i e r a t t h e d r a i n , d u e t o t h e u s e o f a t h i n n e r i n s u l a t o r , a l l o w s h o l e s t o t u n n e l i n t o t h e n a n o t u b e a t h i g h , n e g a t i v e g a t e b i a s . T h e v a l e n c e - b a n d p r o f i l e s f o r t h e t w o ( 1 0 , 0 ) c a s e s , d i s p l a y e d i n F i g . 4 . 1 8 , s h o w t h a t t h i s d o e s n o t h a p p e n w i t h a t h i c k e r i n s u l a t o r . T h i s r e s t r i c t i o n o f t h e o p e r a t i n g r a n g e f o r t h e a g g r e s s i v e l y s c a l e d d e v i c e w a s n o t a p p a r e n t i n e a r l i e r w o r k [109] , w h e r e i t w a s a s s u m e d t h a t a z e r o - h e i g h t S c h o t t k y b a r r i e r c o u l d b e o b t a i n e d w i t h a s m a l l - d i a m e t e r n a n o t u b e . W h i l e t h e i n t r i n s i c s w i t c h i n g t i m e m a y b e a u s e f u l m e t r i c i n t h e e a r l y s t a g e s o f d e v i c e d e s i g n [108] , i t i s w e l l u n d e r s t o o d t h a t p a r a s i t i c c a p a c i t a n c e s b e c o m e m o r e i m p o r t a n t a s d e v i c e s a r e a g g r e s s i v e l y s c a l e d . I n C N F E T s , t h e c o n t a c t s a r e i n e v i t a b l y c l o s e t o g e t h e r , a n d s o i t m a y a l s o b e u s e f u l t o c o n s i d e r a n e x t r i n s i c s w i t c h i n g t i m e w h i c h a c c o u n t s f o r t h e i n t e r - e l e c t r o d e c a p a c i t a n c e s w h e n c o m p u t i n g CG-T h i s i s d o n e b y c o m p u t i n g t h e c h a r g e o n t h e g a t e i t s e l f , r a t h e r t h a n u s i n g t h e c h a r g e o n t h e C N . T h i s i s s t r a i g h t f o r w a r d l y d o n e f r o m t h e f i n i t e e l e m e n t P o i s s o n s o l u t i o n [55] , a n d t h e r e s u l t a n t s w i t c h i n g t i m e s a n d O N / O F F r a t i o s a r e s h o w n i n F i g . 4 . 1 9 . E v e n t h o u g h m o d e s t d i m e n s i o n s h a v e b e e n 4 7 101 10° Q. c 10"1 '"io° 101 102 1 03 1 04 1 05 ^ O l / ^ O F F Figure 4.16: Intrinsic switching time, and ON/OFF ratio for three model CNFETs [55]. employed for the source and drain metallization, which would probably result in excessive series resistance, the contributions of the gate-source and gate-drain capacitances to CG are huge, and lead to much larger switching times than predicted for the intrinsic cases. Numerically, average values of 0.2 and 6 aF are obtained for the intrinsic and extrinsic gate capacitances, respectively. These results are placed in the context of modern devices by comparing them with calculations of the extrinsic switching time of present-day silicon MOSFETs. SPICE simulations1 were performed for a minimum-size NMOSFET using parameters for a commercially available2 90 nm technology. An NMOS silicon device was chosen as its higher current drive makes it superior to PMOS devices [108]. In this way, a fair comparison can be made with the CNFETs studied in the present work, which employ Pd end contacts, which are known to give low barrier heights and, consequently, higher ON currents than are presently possible with n-type CNFETs. Simulations were performed at VDS = IV, for which this technology is optimized, over a VGS range of —1 to +2 V, from which the I O N / Z O F F ratio was obtained for VDD = IV. The total gate capacitance, as returned by SPICE as the capacitance CGG at the chosen operating point, was used to estimate the extrinsic switching time. The result is shown by the dotted line in Fig. 4.19. It is clear that this device out-performs the small-'SPICE simulations were performed by D. Ho of the University of British Columbia. 2STMicroelectronics, http://www.st.com/stonline/. 48 3.5 3 2.5 o 2 x 10 . . . (22,0), T =8nm OX (10,0), T =8nm ' v ' o x (10,0), To=2nm \ 'v. -0.8 -0.6 -0.4 -0.2 Figure 4.17: CN charge as a function of the gate voltage for three model CNFETs [55]. diameter, thin-oxide CNFET that appears to be the most practically relevant of the nanotube devices considered here. Typically, the total gate capacitance of the Si FET is around 200 aF, but superior switching performance is achieved by virtue of the higher drive current and the absence of ambipolar conduction. This suggests that a more extensive study of the switching characteristics of CNFETs, involving different geometries and doped-contact, rather than Schottky-contact, devices needs to be performed before any claims can be made as to the superiority of CNFETs over silicon MOSFETs. This result is in contrast to that of Guo et al., where the sub-picosecond intrinsic switching times of CNFETs were compared to a silicon device [109] using an effective oxide thickness measurement to calculate the gate capacitance [113]. Their conclusion was that CNFETs would actually out-perform the silicon device [109]. 49 0.5h | > 5 I SI _, , , , l _ - 5 0 5 10 15 20 25 Distance from source (nm) Figure 4.18: Valence bands for two oxide thicknesses showing the thinner tunneling barrier for holes in the thin-oxide case [55]. 10 10" 10' fi 10 [ 10 10 10 10 10 . . . . (22,0), Tm=8nm . . . (10,0), 7"M=8nm ox (10,0), T =2nm v ' OX 90nm Si NMOSFET 10' Wo 10° 10* 10° Figure 4.19: Extrinsic switching time, and ON/OFF ratio for three model CNFETs, and a present-day silicon MOSFET [55]. 50 Chapter 5 T h e T i g h t - B i n d i n g M e t h o d A limitation of the effective-mass approximation is its reliance on sufficiently low bias, such that the parabolic approximation to the dispersion relation holds. In order to improve on this model, consider a tight-binding Hamiltonian, and the non-equilibrium Green's function (NEGF) formalism [94,114]. This method allows for non-parabolicity in the bands, as shown in Fig. 4.2, and continues to make the independent-electron approximation. In addition, all of the energy bands may now be considered simultaneously. This will allow for the inclusion of inter-band tunneling without resorting to an energy-dependent effective mass [63, 76], such as in the Flietner dispersion relation [115] with a complex band structure. For quantum wells, introducing an energy dependence into meff may lead to erroneous results [116], so that method should be used with caution. While the NEGF approach has previously been employed to model CNFETs [56-58,65], here the azimuthal potential variations are explicitly included, in contrast to Ref. [56], and a quadratic matrix equation is solved, in order to obtain the required self-energy matrices [117], instead of the usual recursive algorithm [118] which has been used in Refs. [56-58]. In addition, an extension is made to previous works by deriving an expression for the inter-atomic current, which may have use, for example, when considering structures with defects. The key to this method is the computation of the Green's function for a tight-binding Hamilto-nian, so first, a discussion regarding the relevance of this function in terms of computing the charge on the CN surface is presented. 5.1 Relevance of the Green's Function In this section, the derivation presented in Ref. [119] is followed in order to derive the local density of states (LDOS), i.e., the DOS as a function of position, for some atomic system via the Green's function. 51 Suppose that a solution to the time-independent Schrodinger equation, \Z - H(r)]9(r) = 0, (5.1) is required, where Z is a complex number. For some particular choice of Hamiltonian and boundary conditions, this can, in principle, be solved for the wavefunction. One technique for solving differential equations is through the use of a Green's function, G [120]. The appropriate equation to solve is The Green's function, then, gives the response at some point r due to an excitation at the point r', and this fact can be used to solve for $(r) in the entire domain, subject to the boundary conditions and any internal forcing functions. G(r,r';Z), however, also contains information that, when used with the NEGF formalism, obviates the need for an explicit calculation of ^(r). Let &Cti(r) be the z t h eigenfunction of H(r). Using completeness, G may be written as plus an integral over the eigenfunctions in the continuous spectrum [119]. Eq. (5.2) becomes \Z-H{v)]G{v,v';Z)=6{v-v'). (5.2) G = 53/Ci*Cli(r), where Ci are constants. Here, implies a sum over the eigenfunctions in the discrete spectrum (5.3) Now, multiply both sides of Eq. (5.3) by and integrate over position. If Ei is the eigenvalue corresponding to \I/C,i, the result is, after using the orthonormality of the eigenfunctions, *iii(r>) = Ci(Z-Ei). Therefore, the Green's function may be written as (5.4) Now, define (r,r';£,) = limG(r,r ';J5,±iO, (5.5) 52 which is useful since there are poles in Eq. (5.4) at Z = Ei, where Ei is real, and it is convenient to shift them off of the real axis [119,121]. Further, let GQ(r,r';E) = G+(r,r';E)-G-(v,r';E). (5.6) According to Ref. [121], on page 426, lim —]— = p.v.— =F i7r<$(f), (5.7) £->o u ± i£ v where p.v. denotes the principal value. Strictly speaking, this "identity" is an abbreviation for an integral identity resulting from complex analysis, and involves an explicit test function. Here, the full derivation is neglected, and Eq. (5.7) is treated as a true identity. Using Eqs. (5.4)-(5.7) yields G0(r, r';E) = -2m £ '6{E - E^^{v)^'), where the delta function is either of the Kronecker or the Dirac variety, depending on whether the discrete or continuous spectrum is being considered, respectively. Defining D{r;E) = -±:G0(r,r;E), (5.8) it may be observed that this will represent the density of states at r, and if an integration of D(r; E) over all r is performed, the total number of states at a given energy emerges. Note that, if spin degeneracy is included, D should be multiplied by two. The density of states allows for the computation of the carrier distributions in the NEGF formalism, and hence, the wavefunctions do not require explicit calculation. Note that, from Eqs. (5.4) and (5.5), GHv,r';Z) = G(v',v;Z% G-(r,r';E)=[G+(r',r;Erf. As a result, Eq. (5.8) can also be written as D (r ;£) = T - I m { G ± ( r , r ; £ ) } . A generalized LDOS is defined as i \G+ - (G+) +l A D(ry;E).= - L _ ^ _ . U _ (5.9). 53 where A is known as the spectral function. If this function is discretized into matrix form, the equilibrium electron concentration, including spin degeneracy, may be written as n 2 D = -diag /oo Af(E--oo •EF)dE (5.10) where the diag(-) operator extracts the diagonal elements, i.e., the elements corresponding to r = r', of its matrix argument. Out of equilibrium, states that are coupled to a particular lead are filled according to the Fermi function in that lead [114]. Here, a "lead" is defined as a region with constant potential energy in the transport direction, as opposed to the "scattering region" where arbitrary band bending is permitted. As shown in Refs. [94] and [114], a scattering region connected to two semi-infinite leads has a Green's function given by G+ = 1 lim [(E + iO / - HS - E L - En]" 1 , where I is the identity matrix, the S subscript denotes that the Green's function provides information about the scattering region only, Hs is the Hamiltonian for the scattering region, E L and T,R are called the self-energy matrices for the left and right leads, respectively, and A is the atomic density given by 4 3a2v/3 This A factor comes about by assuming that the discrete wavefunction centred on an atom is constant over the area A, and is used in order to give Gj the correct units. Unless otherwise specified, a set of units is employed such that A = 1. As a result, [Gt]'1 - [(GsV'f = 4 - E L + 4 - E f l > (5.11) where £ has gone to zero when we took the inverse. Note that Hs = Hg, and the fact that the Hermitian conjugate of an inverse is equal to the inverse of a Hermitian conjugate, have been used. Eq. (5.11) is often written [^-[(Girf = i ( r L + r f l ) , where Ti = i(Ej — EJ) is known as a broadening function [114]. Performing a simple matrix multi-plication allows Eq. (5.9) to be written as A - J _ 2TT ~ 2?r G j r L (G+)f + G+r* (G+)f] = i - [AL + AR], 54 where Ai = GgVi (Gj)* is the spectral function resulting from the coupling to the lead i. Out of equilibrium, Eq. (5.10) becomes n 2 D -diag / ALf(E-EFL)dE + ARf{E-EFR)AE , (5.12) \J—oo J—oo where EFL and E F R are the Fermi levels in the left and right leads, respectively. It now remains to calculate the Green's function and the self-energy matrices. 5.2 The Hamiltonian Of course, if the Green's function is to be calculated, a suitable Hamiltonian must be defined. For infinite, periodic structures, it is convenient to employ Bloch's theorem and the tight-binding approximation [66], as discussed in Chapter 2. For a finite or semi-infinite crystal, where periodicity breaks near the ends, the similar Hiickel approach [122] may be used. This could also apply to periodic crystals with aperiodic potentials applied to them. In such aperiodic structures, each atomic site should be considered without appealing to Bloch's theorem. Of course, in the semi-infinite case, this results in a matrix of infinite size. However, as observed earlier, if the Green's function can be computed for this system, then the LDOS has also been obtained. Since a finite matrix may be inverted directly, this case will not be considered. In a discretized form, the equation for the Green's function, without appealing to self-energies, is given by [94] G± = \\m[{E±\i)I - H)~l . The Hiickel Hamiltonian looks very much like the tight-binding Hamiltonian derived earlier [122]. Using a suitable energy reference, considering a (2,0) CN for specificity, the Hamiltonian may be 55 written as H = 7 7 U7 7 7 7 7 We 7 7 (5.13) t/4 7 7 U3 7 7 7 ^2 7 7 7 t/i I such that the end of the tube is represented by elements in the bottom-right corner of the matrix. Note that the diagonal elements, Ui, are given by the local potential energy. Of course, the pattern that has been shown .here is replicated along the diagonal out to infinity. The exact form of the matrix depends on the ordering of points, but this is the result for one particular ordering. Since the Hiickel derivation is similar to the tight-binding case, it is simply noted that by explicitly accounting for each atom, a matrix of a similar form as in tight-binding, but with no wavevector dependence, is obtained. The price is that there is now an infinite matrix to consider. Details may be found starting on page 627 in Ref. [122]. 5.3 Computing the Green's Function In this section, a technique for computing the Green's function [117] is presented. For an infinitely-long CN, a general version of the nearest-neighbour, tight-binding Hamiltonian of Eq. (5.13) may be written in the block form H r t T g r where it may be recalled that a nanotube has periodicity in its atomic positions. Indeed, the following method could, in principle, be used for any crystalline material. The periodicity is defined by the 56 CN primitive unit cell, as discussed in Chapter 2, containing 2N atoms. Here, r represents the coupling between these primitive unit cells, T represents both the local potential energy and the coupling between atoms within a primitive unit cell, and the subscripts L , S, and R represent the left lead, the scattering region, and the right lead, respectively. Note that all of these quantities are matrices, and that the left and right leads extend off to infinity. The magnitude of the atom-atom coupling is given by |7|, and in this case, a nominal value of 2.8 eV is chosen as appropriate for the pz orbitals [67]. In addition, note that rows corresponding to each unit cell of the leads repeat along the diagonal of H out to infinity due to the form of the potential energy there. Only potential energy variation normal to the transport direction is permitted in the leads, and any arbitrary potential variation is permitted in the scattering region. The Green's function, which contains information about the LDOS [94,114], is given by G=(EI-H)-\ where the notational change G = G+ has been used, and the appropriate limit has been implicitly taken. Since the scattering region is the only relevant part of the calculation, the infinite matrix may be truncated using the method described in Ref. [94]. This allows the Green's function for the scattering region to be computed via G S = ( ^ / - T s - E L - E f l ) - 1 , where E L = T + G L T , • (5.14) ^R = TGRT^ (5.15) are the self-energies, and GL and GR are the Green's functions for the left and right leads, respec-tively. The self-energy matrices are finite, and as a result, the inverse may be computed directly, or by using a recursive algorithm [118]. For the self-energy calculations, two coupled, linear matrix equations may be solved simultaneously, as in Ref. [56]; however, a more efficient method is to solve a single, non-linear matrix equation. Due to the tight-binding form of the Hamiltonian, the self-energies may be written in block form 57 as E r = C L 0 0 b o o o o o o o o Efl= 0 0 0 0 0 where the diagonal elements, from upper-left to lower-right, have sizes 27V x 2N, (ns — 47V) x (ns — 47V), and 27V x 27V, respectively, where ns is the total number of atoms in the scattering region, and it has been assumed, without loss of generality, that ns is an integer multiple of 27V. The LDOS is given by the diagonal elements [94,114] of GS ( S i - E i + E a - E ^ ) Gs, and if 67s is written in block form as Gn G\2 G13 Gs = G21 G22 G23 G31 G32 G33 Gn G\2 G\z I 0 0 G21 G22 G23 = 0 0 0 G31 G32 G33 0 0 I it may be observed that the second block column of Gs is irrelevant. As a result, it is sufficient to solve (EI — Ts — E L — E F I ) for Gij. The LDOS may, then, be computed as D(r; E) = diag + 4)6^ + G i 3(te + 4)Gj3] , where i = 1, 2, and 3 for atoms in the unit cell adjacent to the left lead, away from either lead, and adjacent to the right lead, respectively. It now remains to describe how one may compute G L and GR in order to obtain the self-energies from Eqs. (5.14) and (5.15). 58 5.3.1 G r e e n ' s F u n c t i o n f o r a S e m i - I n f i n i t e L e a d Previously, Appelbaum et al. have presented a method for calculating the required elements of GL and GR for the case where r is a scalar [123], and Venugopal et al. have presented that for the case where r = rt [124]. In the former work, the recursive algorithm [118] is implicitly employed to arrive at the solution, and this recursion may also be used for the general case considered here, as in Ref. [56]; however, as mentioned previously, by framing the problem as a quadratic matrix equation, the computation may be performed more efficiently. GL is the inverse of E I - H L = •rt EI - TL T •ft EI - T z T EI — TL Due to the infinite size of this matrix, this equation can also be written compactly as E I - H L = GL is defined as r m gu 9si 9s where gs are the surface elements of GL-Eq. (5.16) may also be written as EI-HL -rt EI — TL G L = EI-HL T rt EI - T L (5.16) EI-HL T 9i 9is I 0 rt EI - TL 9 si 9s 0 I This matrix equation can be solved to get gs = [EI - T L - T\EI - HLrlr]_1. Of course, this still requires the inversion of an infinite matrix. However, the nearest-neighbour, tight-binding structure of r may be exploited in order to get that T^(EI — HL)~1T = T^gsr. As a result, the problem now is to solve the nonlinear equation gs = [EI-TL-^gsT]-\ 59 (5.17) which may be most simply done using a fixed-point iterative scheme, as in Ref. [123], such as ^+i) = [ E / _ T z _ _ r t 5 ( i ) r ] - i j where the superscript indicates the iteration number. This is a form of the recursive algorithm for self-energy calculations [118]. While this method does appear to converge for all cases attempted thus far, the convergence of the fixed-point, or recursive, method can be quite poor. As a result, the use of Newton's method on Eq. (5.17) is considered. The residual, r, is defined as r{ga) = T^gsTgs - (EI - TL)gs +1, which is equal to zero when ga is the exact solution. This is a quadratic matrix equation, and the error matrix, e, that may be applied to some approximate solution, gs, is sought, i.e., r(gs +e)= r\ga + e)r(g3 + s) - (EI - TL)(gs + e) + I is solved for e, where r(gs + e) = 0. This leads to r(gs + e) = r(g3) + FSa (e) + Sere, where Fg, (e) = r h r g s + [ T ^ . T - (EI - TL)]e (5.18) is related to the Frechet derivative of r(ga) in the direction e. The. Newton method is defined by solving r(gs) + Fgs(e) = 0 for e. As noted in Ref. [125], this is a special case of the generalized Sylvester equation, and may be solved using the method described in Ref. [126]. Due to the sensitivity of Newton's method to the initial guess, iterations are damped using exact line searches [125], where, instead of performing the usual Newton update of g s l + 1 ^ = gi^ + £, the modified form of <7^ +1) = + -de is used, where •& is a parameter chosen such that the square of the Frobenius norm of r(g^ + i9e) is minimized at each step. Details of this procedure are provided in Ref. [125], and in Appendix E. In general, exact line searches are required in order to obtain convergence with Newton's method. Due to the increased computational effort per iteration for this method, the iterative solution begins, in practice, with fixed-point iterations for a while, and then switches to Newton with exact line searches. Typically, only a few Newton iterations are required in this hybrid approach. 60 A completely analogous procedure holds for the calculation of the surface elements of GR, except that r is exchanged with T*, and vice versa, in all of the above derivation. 5.3.2 Test Cases The simplest case with which to check this method is for the case where the local potential is some fixed constant. This case is a periodic structure extending out to infinity in both directions, and therefore, the Bloch condition may be applied in order to compute the DOS using a single, primitive CN unit cell. Equivalently, one may employ the zone-folding scheme in order to obtain the DOS from the dispersion relation for graphene (see Refs. [4,6], for example). Both of these methods are described in detail in Chapter 2. In this case, TL = Ts = TR, and the NEGF method is used to compute the LDOS from Gs [94,114]. In order to numerically handle the van Hove singularities in the CN DOS, energies were perturbed by a small, positive imaginary number, i£. For the results shown here, the perturbation was £ = 2 x 10~6, and convergence was declared when the infinity norm of r(gs) became less than 10 - 6. As shown in Fig. 5.1, this method produces the correct form. •ni • r b » > « r i — = • 1 - 1 . 5 - 1 - 0 . 5 0 0 . 5 1 1.5 E n e r g y ( e V ) Figure 5.1: Comparison of the standard tight-binding density of states (solid line) with that com-puted using the self-energy method (o) for a (16,0) carbon nanotube [117]. 61 The second case to consider is that of a finite barrier. In this case, quantum mechanical reflection for carriers above the barrier, and the possibility of tunneling through the barrier, are expected. These phenomena should modify the form of the LDOS at each position. Fig. 5.2(a) shows the LDOS for this case, and the solid lines show the conduction and valence band edges shifted by the local potential energy. Regions of high electron concentration are shown as bright patches, and since the Green's function represents outgoing waves at infinity, a vanishing density is obtained in the scattering region for energies lined up with the bandgap in the leads. Also note that the phenomena of inter-band tunneling, which may have important technological applications [63,127], is clearly shown in the figure. The jumps in density with energy are due to the contributions of higher bands, and intra-band tunneling is also evinced by the penetration of carriers into the bandgap region. It is noted, here, that in Ref. [114] a subsequent summation was required for each transverse mode, or band, and that these modes were assumed to be independent. In this case, a single calculation provides a complete multi-band result with allowed inter-band transitions. Finally, the LDOS for a potential that varies both azimuthally, and along the length of the CN, is computed. The longitudinal potential variation is multiplied by 1 + O.2sin(0). Fig. 5.2(b)-(d) shows the resultant LDOS at <fr = 0, TT/2, and 37r/2 radians, respectively. While the variation with angle is readily identifiable, particularly in the valence band, it is interesting to note that the local band structure appears to be more influenced by the average potential around the circumference than by the local potential. This may be seen by the superimposed conduction and valence band edges, which have simply been shifted by the local potential energy. The LDOS itself does not undergo such a shift directly; rather, its magnitude is modified, with the gap apparently taking the position appropriate for the average energy shift, i.e., the shift at cj> = 0 in this case. This may bode well for methods that utilize an average potential on asymmetric structures [55,56]. Careful examination reveals that the mid-length conduction band LDOS is highest for r/> = 37r/2, which corresponds to the angle with the lowest potential energy, followed by 0 = 0, and then <p = TT/2. This is the expected result; the unexpected result is that the LDOS is not rigidly shifted by the local potential energy. While this method yields excellent results, further analysis is required due to the possibility of a multiplicity of solutions that may satisfy the quadratic matrix equation. It is well-known that equations of this type may have any number of solutions, including no solution at all. At present, it can only be stated that there is at least one solution to the equation, and that it is a physically 62 10 15 20 Distance (nm) 25 30 1 0 15 20 Distance (nm) 25 30 (a) (b) Distance (nm) Distance (nm) (c) (d) Figure 5.2: Local density of states for longitudinal and azimuthal potential variation. Shown are: (a) a barrier potential with a maximum shift of 0.8 eV longitudinally, and no azimuthal variation; the same barrier potential as in (a) multiplied by 1 +O.2sin(0) and evaluated at angular positions of (b) 0, (c) 7r/2, and (d) 37r/2 radians. Solid lines show the result of shifting the conduction and valence bands by the local potential energy. 63 reasonable one. 5.4 Current Operators Ultimately, the current-voltage relationship is of interest. Consider, then, how the real-space current operator may be derived when considering quantum transport. Both effective-mass and tight-binding Hamiltonians are considered in order to better justify the tight-binding case by comparing it to a well-known result. One- and two-dimensional operators will be derived on a Cartesian grid, and the conclusion will involve the derivation of an expression for the current in a particular CN with a tight-binding Hamiltonian. The CN case is treated separately since the Hamiltonian allows for coupling between atomic sites that are not arranged on an orthogonal mesh. As a result, the CN matrix will differ from that of atoms on a rectangular or square crystal lattice. The derivation begins with the time-dependent Schrodinger equation, <9\I> ih— = H9, (5.19) where t is time. In addition, there is the conjugate equation, i ^ = -{H9)l (5.20) Following the approach of Ref. [128], note that dtK ' dt dt Using Eqs. (5.19) and (5.20), ~(^^>) = ^ [ V # # - . (5.21) An expression similar to this occurs in hydrodynamics if the associations V = V- J = - [H^V] , (5.22) are made, where V is like a fluid density, and J is a probability current density. The tilde is a reminder that these quantities are functions of energy. In order to obtain the total probability current density J , J must be integrated over E. To convert J to an electronic current density, simply multiply by the charge of a carrier q. It now remains to discover the appropriate expression for J . The goal is to write the right-hand side of Eq. (5.22) as a divergence expression. 64 5.4.1 Effective-Mass In the effective-mass approximation, the explicit form for the Hamiltonian is H = --^-V2 + U(r). 2meS Since the potential energy is real, Eq. (5.22) can be written as V • J = —— [*tv2* - (v2*t) 9] . Z17Tl e ff If it is recalled that V 2 * = V • (V*), this becomes V • J = —^— f* fV • (V*) - V • (V**)*] . 2imeff A standard vector identity is V • ( ^ V * ) = tf'V • (V*) + V * • Vtft, and similarly, V • (tfVtft) = * V • (V* f) + • V * . Subtracting these identities implies that v J = -T^—V • IVv * - w i r /T l . 2irneff L 1 Therefore, the current density can be defined as h 2im, k [* +V* - * W f ] = Re •eff WTlefF $ty$ (5.23) Note that this expression is not unique. One could, for example, add any function to J that has a vanishing divergence, and still obtain Eq. (5.22). Here, an appeal is made to the correspondence principle, which states, from page 29 of Ref. [128], "Quantum Theory must approach Classical Theory asymptotically in the limit of large quantum numbers." In order to satisfy this condition, "one establishes in principle that there exists a formal analogy between Quantum Theory and Classical Theory." In this case, it is noted that the formal analogy of the classical momentum is the operator (h/i) V. As a result, for a plane-wave wavefunction, the form of Eq. (5.23) is analogous to the velocity multiplied by the density, as in the classical current. An additive term with a vanishing divergence is, therefore, not necessary. 65 5.4.2 Tight-Binding In the tight-binding approach, the right-hand side of Eq. (5.22) will, again, be written as a divergence expression. As in the effective-mass case, it will be assumed that the current density is given by this expression without adding any function with a vanishing divergence to the result. Previously, an appeal was made to the correspondence principle. Here, the appeal is made to the discrete version of Eq. (5.23), and no attempt is made at trying to present a formal proof. It is assumed, throughout, that the atomic spacing, a, is uniform. The first step is to discretize Eq. (5.21). If #3 (*l = [*j * 2 * J • • • * j v ] , where * i is the wavefunction at the ith atom, and N is the total number of atoms in the system, then Eq. (5.19), and its dual equation, may be written as . i f i^ |*)=H|*>, (*| = -<*| Hi l*> = Proceeding as before, | (I*) <*l) = ^ (H\9) <*| - |*) <*| ^ ) , (5.24) where it is noted that the diagonal elements of this matrix equation are the discrete equivalent of evaluating Eq. (5.21) at the lattice sites. In order to make the right-hand side of Eq. (5.24) look like a divergence expression, it is helpful to recall the formal definition of the divergence [129], lim A - . 0 I ^ J . n d " Y , (5.25) where m is one less than the number of spatial dimensions of J , A is the volume enclosed by the surface s which surrounds the point at which you are calculating the divergence, and fi is a unit 66 vector that points in the outward-normal direction to the infinitesimal area element dmr. If working in dimensions different than three, the appropriate generalizations of the terms volume, area, and surface should be employed. Cartesian Lat t ice First, consider the case where the lattice of atoms forms a Cartesian grid, as in Fig. 5.3. The tight-binding Hamiltonian is assumed to include only the effects of the nearest neighbours. • • • • • • • • • • • • • • • • • • • • Figure 5.3: Sample one- and two-dimensional lattices on the left and right, respectively. In one dimension, the Hamiltonian can be written as a tri-diagonal matrix of the form 0 7 0 ••• 0 0 0 7 0 7 0 0 0 I H I 0 0 0 7 0 7 0 0 0 ••• 0 7 0 where all the real, diagonal elements; such as the local potential energy, have been eliminated since they will cancel out in the diagonal elements of Eq. (5.24). Substituting the ID, tight-binding Hamiltonian into Eq. (5.24) yields diagonal elements, away from the ends, of the form A ID approximation to Eq. (5.25) is 67 are defined, where J i + i is the magnitude of the current density between atoms i and i + 1. It now remains to compare these expressions with the effective-mass result. Using the finite-difference approximation, one can view the tight-binding Hamiltonian as an approximation to the effective-mass Hamiltonian with h? 7 = - T r — • (5-26) 2a2meff Using finite differences on the effective-mass current expression leads to 1 h Ji+h 2irneff h 2 2iames where it has been assumed that the wavefunction between sites i and i + 1 is equal to the average of the wavefunctions at the nearest atomic sites. Using Eq. (5.26), it becomes apparent that the tight-binding result is equivalent to the effective-mass result. Of course, this equivalence is expected since the tight-binding Hamiltonian may be viewed as an approximation to the effective-mass Hamiltonian when working with atoms on a uniform, Cartesian grid. Using a 2D, tight-binding Hamiltonian, the result is that where a double subscript to the wavefunction, ^ fij, is now required in order to denote the ith site in the x-direction, and the j t h in the y-direction. The 2D approximation to Eq. (5.25) is V J where the surface, s, is defined as a square with area a2 centred on the site at Note that it is assumed that the normal component of the current density is constant along each of the four sides of s, so the integral is just the normal component of the current density multiplied by a. This leads to similar definitions for J, 68 where Ji+1 ^ is the magnitude of the current density component in the y-direction crossing the plane between rows i and i + 1, and in the j t h column, of atoms. J J J + I is defined similarly. By direct analogy to the I D case, it is concluded that this expression is also equivalent to the effective-mass case. Hexagonal Lat t ice In a hexagonal lattice, such as in a CN, there is a two-point basis. As a result, two different equations emerge, depending on which atom in the basis is being considered. The primitive unit cell of the hexagonal lattice is a parallelogram [66]. If this unit cell is divided into two equilateral triangles, a convenient s over which to approximate Eq. (5.25) is found. The hexagonal lattice, primitive unit cell, and s are shown in Fig. 5.4. Figure 5.4: Hexagonal lattice (solid) with a primitive unit cell (dashed parallelogram) divided into two equilateral triangles s. The appropriate approximation to Eq. (5.25), using the same assumptions as before, and con-69 sidering the higher of the two atoms, denoted by a superscript h, in the primitive unit cell, is 4 V 3h ~ — (Jt + J\ + Jy) where the subscript indicates the direction of the surface through which the current passes. Similarly, for the lower atom, denoted by a superscript I, 4 Utilizing the tight-binding Hamiltonian, dt L K 1 1 ift Therefore, A [<p'(«.')tl = 1 r*|_(*')t + ($')t + ' (*')t_ a t ift 4ift is defined, and similarly for the other current densities. 5.4.3 C u r r e n t O p e r a t o r S u m m a r y After all of this work, a summary of the current operators, J, is presented, where it is noted that it has been implicitly assumed that <$> has been normalized to give the electron density at a given atomic site, i.e., if i/j^tp gives the electron concentration per atom, then fy^ty = ip^ip/A. In this section, the A normalization is explicitly included. 1. The effective-mass current density is given by J = J^L (tfVtf* - * f V * ) . 2. The I D , tight-binding current density operator is given by J^-'-A1 ' i = < + 1 > (5-27) l f t 1 0 , + 70 where Jij is the element in the z t h row and written as J = diag(j|$) where J = 3. In two dimensions, assume that the elements !#> = j t h column of J. The current density may be + (5.28) J 3 2 2 Ji 2 of \9) are ordered by increasing x first, i.e., # 2 , 1 # 1 , 2 # 2 , 2 where 7V X is the number of points in the ^-direction. With this ordering, the current operator for the 2>direction is also given by Eqs. (5.27) and (5.28), where J is ordered similarly to |#). The only difference is that, for a finite number of atoms in x, the matrix will have only zeros in the rows corresponding to the last atom in that direction. A similar expression holds for the y-component of J . 4. The current operator for a tight-binding Hamiltonian on a CN hexagonal lattice depends on the assumed chirality due to the periodicity in the circumferential direction, and it also depends on which atom is being considered in the two-point basis. Rather than trying to present all of the various permutations that are possible for a general CN, a demonstration of how one may compute the current through the cross-section of the CN, depicted in Fig. 5.5, is presented. From the current density expression, it may be seen that the current across the dashed line will be related to the wavefunctions at only four atomic sites. If a subscript denotes which 71 Figure 5.5: Atomic structure of a (2,0) carbon nanotube for which a calculation of the current across the dashed line is desired. Note that the atoms on the top of the figure are identical to those on the bottom due to the periodicity of the tube. atomic site is being considered, the current, T, may be written as where spin degeneracy has been accounted for in the normalization of the wavefunction, and the tilde, again, implies that this quantity should be integrated over energy to give the total current. In order to get this result, it has been assumed that the current is constant over each face of the triangle depicted in Fig. 5.4, and an integration has been performed over the dashed line in Fig. 5.5. 5.5 A l t e r n a t e C u r r e n t E x p r e s s i o n In the previous section, current operators were derived that provide the current density spectrum as a function of position and energy. While this will be interesting in terms of being able to inspect the path of carrier flow in the device, the terminal characteristics are of primary importance. In this case, a compact expression may be derived for the terminal current following the approach detailed 72 in Ref. [130]. This compact expression is commonly used in mesoscopic physics. The Schrodinger equation, in matrix form, is {E + i$)I-HL -r 0 -rt (E + %)I - Hs -T #s = 0 -rt (E + i$)I-HR 0 where 9L,S,R = ^L,s,Rl'V~h is the wavefunction in each of the three regions under consideration, i.e., the left lead, the scattering region, and the right lead, respectively; and it is recalled that the Green's function involves the limit as £ approaches 0 +. Using ipi from this point on, it is convenient to let ipL — <PL + XL and ipR = <PR + XR> where <pi satisfies [(E + %)I-Hi]<pi = 0, and represents the wavefunction in the lead before coupling to the scattering region. After coupling, the unperturbed wavefunction is modified by Xi- Performing some algebra yields i>s = \{E + \i)I-Hs-T,L- VR}-1 [TVL + TVR] = Gs [TVL + T<pR] , (5.29) XL = GLTi>s, (5.30) XR = G R T ^ S - (5.31) From Eq. (5.24), trace (^ps^l) = ^ T R A C E ( T V L V ' S - ipsip{T + TipRips - tps^RT^ = - / V • Jd 2r, since Hs = Hs, and the trace, i.e., the sum of the diagonal elements, is not dependent on the order, of the matrix multiplications in its argument. Note that taking the trace over the scattering region elements is the discrete equivalent of integrating over the scattering region S up to a normalization factor of A. This factor cancels out when working with ipi instead of Utilizing the divergence theorem, and noting that current may only travel out through the leads, yields -1L +IR = trace (T^L^S ~ i^s^r + TipRips ~ V'sV'Ji'rt) . where 11 and TR are the currents out of the scattering region at the left and right interfaces. Clearly, TL = TR in steady-state, since the time derivative of ip^ip vanishes. Of course, this implicitly assumes 73 that all of the charge in the scattering region is due to the leads, and hence there is no recombination or generation allowed in the scattering region. Considering only the left lead, XL = T|trace ( T V L V ' S + T^XL^S - ^s9LT ~ i>sXLT) • In order to write this in a more convenient form, substitute Eqs. (5.29)-(5.31) into this expression. If it is noted that [130] t fiAiX where fi is the distribution function in the ith lead, the result is, after some algebra, XL = |trace {fLTLAR - fRYRAL). This may be simplified even further if it is noted that [130] trace (TLAR) = trace (TLGS^RG^J = trace (GsTLGlrR) = trace (TRAL). The final expression for the current is [114,130] XL = ^trace [rLAR (fL - fR)], (5.32) which should be multiplied by two if spin degeneracy is included. The same expression results from considering the right contact. In this form, Eq. (5.32) bears a striking resemblance to the Landauer current formula [93]. If the transmission probability for the scattering region is associated with trace (TLAR) [114,130], then the formulae are identical. Use of this compact expression is convenient since only diagonal elements of matrices that are already known from the charge calculation are required. In order to calculate the inter-atomic current density, further computation of off-diagonal elements would be needed. 5.6 S o l u t i o n P r o c e d u r e In order to arrive at converged solutions, two algorithms have been implemented. The first is related to the standard Gummel technique as described, for example, in Ref. [131]. This method, when it converges, is useful due to its good convergence rate [56,58]; however, it is found that a good solution 74 is actually less likely to be found with this method than with an alternate Newton iterative scheme with an approximate Jacobian [132]. Doped-contact CNFETs are considered, with a potential energy in the scattering region that joins smoothly to that in the leads. The scattering region, thus, includes some length of the doped contacts. The Gummel iterative scheme proceeds as follows: 1. Ensure that an integral number of primitive unit cells are in the scattering region. This is really just a bookkeeping issue, and ensures that the leads have symmetric structures. This is useful when computing coupling matrices. Since the leads have no potential variation in the transport direction, the amount of doped contact included in the scattering region may be increased without changing the solution. 2. Tetrahedralize the domain for the finite element Poisson solver. In practice, this is done using the freeware program TetGen [102]. 3. Compute the Hamiltonian, and atomic positions, for atoms in the scattering region. The Hamiltonian is a sparse matrix with non-zero entries only for nearest-neighbour atoms. The atomic positions are needed in order to interpolate between the atomic grid, used for the Green's function calculation, and the Poisson grid. 4. Compute the Hamiltonian for one CN unit cell of the leads. This is one block of the matrix computed in the previous step. 5. Compute the stiffness matrix for the finite element Poisson solution. Details of this calculation are provided in the Finite Element Electrostatics section of Chapter 4. 6. Guess initial quasi-Fermi levels for holes and electrons, E p v and E p n , respectively, for use in the nonlinear Poisson equation detailed in Step 7. 7. Solve a simplified nonlinear Poisson equation. This is done to improve the convergence prop-erties of the overall iterative scheme [56,58,131]. The idea is to solve W = -±5(p - Rt) (P2D(V) - n 2 D(V) + Nf), where (5.33) 75 n 2 D = r / ^ d £ V ( 5- 3 4) 7o e X p ( - ^ ) + l and D(E) is the 2D density of states, based on symmetric electron and hole bands. In practice, the analytic D(E) of Ref. [133] divided by the CN circumference is used. This Poisson equation may be straightforwardly solved using Newton's method. Note that, on subsequent iterations, EpP and Epn are calculated such that P2D(10 and n2T>{V) are equal to the carrier concentra-tions computed via the Green's function. If the charge and potential are self-consistent, the solution to this nonlinear equation will converge in the first step. 8. Solve for the Green's function in the scattering region in order to obtain the charge there via Eq. (5.12). Note that the net charge density is given by q(l — n2o) where the first term represents the ionic charge left over if no pz electrons are present. In order to facilitate calculation of the integrals in Eq. (5.12), the first term can be written as —diag 7T /EH PEN roo ALdE- ALf(EFL-E)dE + ALf{E-EFL)dE •oo J— oo JEN . n 2 D = -IdiagJ£ r WK*-*»W>AL ^ DE and similarly for the second term, where E^ is the "charge neutrality level" defined such that the integral of the total density of states up to this energy yields the atomic density A, where it is noted that the total LDOS is related to AL + AR, and that A = 1. Using this yields 8ga(E-EN(z))AL / - o c exp [sgn(£ - EN{z)) (^f^)] + 1 + r ^{E-EN{Z))AR I > ( 5 3 5 ) J-oo exp [sgn (E - EN{z)) (^f^)J + 1 J where E^{z) is shifted by the local potential energy, and sgn(-) is the sign function. Using the charge neutrality level, and the sign function to obtain more compact and computationally tractable carrier integrals has been previously seen in Refs. [23,56,134], for example. It is sometimes convenient to associate E — E^(z) < 0 with p2o and E — EN(Z) > 0 with n2D, particularly in the Gummel method where quasi-Fermi levels for both holes and electrons may require computation. 9. Compute the quasi-Fermi levels corresponding to the electron and hole concentrations from Step 8. This is done using Newton's method to compute Epp and EFn from Eqs. (5.33) and (5.34) given known p2D and n 2D-76 10. Repeat Steps 7-9 until the charge from Step 8 is consistent with the potential, i.e., the sim-plified nonlinear Poisson equation converges in a single Newton step. The Newton algorithm is identical in Steps 1-5. Following these precalculation steps, the algo-rithm proceeds as follows: 6. Guess an initial potential distribution. In practice, the initial guess is generated by performing Steps 6 and 7 from the Gummel algorithm above, with quasi-Fermi levels chosen such that Epp = Epn — s 0 , in the left lead and the left \ of the scattering region —QVDS J in the right lead and the right | of the scattering region qVDS ,, otherwise 2 7. Solve for the Green's function in the scattering region in order to obtain the charge there. This is done via a charge neutrality level as in Step 8 from the previous algorithm. The approximate Jacobian of the Poisson equation, KV = C[V], is also computed at this time, where it is noted, by the square brackets, that C is a functional of V, i.e., it depends on V everywhere on the CN surface, and not just on the value of V at a single point. The elements of the full Jacobian matrix are given by 11 dVj' where the first term is the finite element stiffness matrix, which contains mostly zero elements. As a result, sparse matrix methods may be used to store, and perform operations with, this term. The second term is also largely zero with nonzero elements only for i and j corresponding to points on the CN surface. Unfortunately, this term is still much too computationally intensive to evaluate fully, and is still too large to store if each CN atom is being treated. As a result, it is approximated by taking the derivative with respect to the contact Fermi levels [132]. This results in a matrix with nonzero elements only along the main diagonal. 8. Update the potential as usual in Newton's method with the approximate Jacobian. 9. Repeat Steps 7-8 until the potential converges, i.e., until the charge density and potential are self-consistent. 77 Due to the time-intensive nature of these computations, it has been helpful to use parallel pro-cessing in order to compute the charge density. This was implemented on the energy integrals using the message passage interface [135]. The domain of integration was split into several subdomains, and these smaller integrals could be computed simultaneously. Upon completion, each small integral was sent to the master process for summation in order to complete the computation of the total integral. 5.7 R e s u l t s a n d D i s c u s s i o n This section begins by considering two device structures, depicted in Fig. 5.6: an asymmetric double-gate (DG) device, and a semi-cylindrical-gate (SC) device. In both cases, the source and drain are (a) (b) Figure 5.6: Schematics of the (a) double-gate and (b) semi-cylindrical-gate CNFET geometries. doped nanotubes of the same chirality as the intrinsic channel, and are shown by the different shad-ing in the figure. The DG case is close to structures presently being investigated experimentally, and the SC device may be a realizable tri-gate-like structure. It is assumed that the entire simulation space outside of the CN is filled with a homogeneous dielectric, taken to be Si02 unless otherwise specified. A Dirichlet boundary condition is applied on the gate, and homogeneous Neumann bound-ary conditions are used on the remaining boundaries of the simulation space, including on the outer ends of the source and drain, as in Ref. [57]. 78 The reason for the Neumann conditions on the source and drain stems from the fact that the devices involve doped, semiconducting, ballistic contacts, and as a result, it is no longer appropriate to use Dirichlet boundary conditions, as was done in the metal-contacted, Schottky-barrier case [136]. Rather, recall that some length of the source and drain are included in the simulation space. Neumann conditions are applied at their ends [114,136], such that the source and drain potentials effectively float at the edges of the domain, and the drain-source voltage is included only when evaluating Eq. (5.35), where EFL = 0 is the Fermi level deep in the source, and EFR'= —qVos is the Fermi level deep in the drain. In order for the solution to be valid, it is required that the local potential becomes constant in the longitudinal direction well before the end of the simulation space. To illustrate this, consider the potential in Fig. 5.7. Shown are the results for gate thicknesses of 10, Distance (nm) Distance (nm) (a) (b) Figure 5.7: Local CN potential for (a) the double-gate and (b) the semi-cylindrical-gate devices with VGS = VDS = 0. The doped regions are 20 nm long at either end, and are doped n-type to a density of 0.38 nm~2. The intrinsic region is 7nm long. Shown are curves corresponding to gate thicknesses of 10nm (solid), 2nm (dashed), and 0.1 nm (dash-dotted). Data is taken along the top of the nanotube. 2, and 0.1 nm by the solid, dashed, and dash-dotted lines, respectively. In each case, the potential is forced to become flat by the imposition of the boundary conditions. In both the DG and SC cases, the potential becomes flat well before the boundary for the 2 and 0.1 nm thicknesses, and these 79 solutions may be considered to be valid. In the 10 nm case, particularly for the SC device, this is not the case, and the solution is not valid. For a valid solution, charge neutrality must be achieved within the simulation space. In order to correct this, longer doped regions must be simulated; however, due to the atomistic approach, this incurs a huge computational cost. At equilibrium, insight may also be gained into the azimuthal potential variation that results from these non-cylindrical structures. This is depicted in Fig. 5.8 for the CN surface potential in the middle of the intrinsic region as a function of azimuthal angle, where 0 = 0 corresponds to the top of the CN, i.e., the closest point to the top gate in the DG case. The DG structure exhibits 0.3 0.25 - 0.2 03 . | 0.15 CL 0.1 0.05 ° -3 - 2 -1 0 1 2 3 <t> (rad) Figure 5.8: Local CN potential in the middle of the intrinsic region as a function of azimuthal angle for the double-gate (dashed) and semi-cylindrical-gate (solid) devices. greater azimuthal variation as it is further from the completely symmetric, cylindrical geometry. In order to investigate the impact of the azimuthal variation, and the geometrical differences, results are presented for five structures: a DG device with zero gate thickness, similar to that presented in Ref. [57]; a DG device with 5nm gate thickness; a DG device with 5nm gate thickness where the axial potential is used in the NEGF solution rather than the proper surface potential; an SC device with zero gate thickness; an SC device with 2 nm gate thickness. In all cases, a p-i-n device is considered, which is interesting due to its potential for an exceptional subthreshold slope [63,127]. The doping density is set to 0.1 nm - 2 in the p- and n-regions. 80 These devices have also been studied using the Flietner dispersion relation [63,127], and by Koswatta et al, using a two-band NEGF simulation [65]. These results are an extension in that they use NEGF to consider all of the bands simultaneously, in addition to the consideration of geometrical variations. The simulation space includes 15 nm each of the source and drain doped regions, in addition to the 7 nm intrinsic channel, and a 2 nm oxide thickness is chosen on top of the CN, and 20 nm below in the DG case. In order to better understand the operation of a p-i-n device, consider the local density of states, and the energy spectrum of the electron density for the SC device. These are shown in Fig. 5.9. The bright patches show regions of high density, the solid lines show the field-free conduction- and 0 10 20 30 0 10 20 30 Distance (nm) Distance (nm) (a) (b) Figure 5.9: (a) Local density of states and (b) energy spectrum of the electron density for the semi-cylindrical-gate device. Solid lines correspond to the field-free conduction- and valence-band edges shifted by the local potential energy, and the dashed lines show the source and drain Fermi levels deep in the contacts. The dashed lines terminate at the ends of the doped regions. Data is taken for VGS = -0.057 V and VDS = 0.4 V. valence-band edges shifted by the local potential energy, and the dashed lines indicate the positions of the source and drain Fermi levels deep in their respective contacts. The dashed lines extend throughout the doped regions. In Fig. 5.9(a), the electronic states that are available to be populated by source and drain injection are shown. As seen in Chapter 4, the jumps in density correspond 81 to the multiple bands present in the CN. Fig. 5.9(b) shows how these states are actually populated by source and drain injection. It appears that electrons may tunnel from the valence band to the conduction band near the doped region in the drain. In order to confirm this, the energy spectrum of the current is plotted in Fig. 5.10, which clearly shows electrons tunneling from occupied states in the valence band to empty states in the conduction band. Note, too, that many of the electrons originate in the second valence band, and flow into Distance (nm) Figure 5.10: Energy spectrum of the current for the same device and conditions as in Fig. 5.9. the first conduction band. This result highlights the importance of being able to perform a full multi-band simulation, and would be missed in Ref. [65], for example, where it appears that only the highest valence- and lowest conduction-bands were considered. Turning to the terminal characteristics of the model devices, the ID-VGS data is summarized in Fig. 5.11. For these devices, the subthreshold slope is significantly worse than the thermionic limit of 60mV/dec, with a value around 300 and 700mV/dec in the semi-cylindrical- and double-gate cases, respectively. In order to investigate this poor performance, the energy spectrum of the current may be utilized. Close inspection reveals that a significant tunneling current exists for all energies between the contact Fermi levels with a maximum for energies where the tunneling distance is shortest, as shown in Fig. 5.10 for example. For these model devices, the current is appreciable even in the OFF state, with a value on the order of 0.3 pA. This OFF-current is the root cause of 82 10"5 -0.2 0 0.2 0.4 0.6 0.8 VGS ^ ' Figure 5.11: TD-VGS for double-gate and semi-cylindrical-gate devices. The D G devices are zero gate thickness (•), 5nm gate thickness (o), and 5nm gate thickness where the axial potential is used (A) . The SC devices are zero gate thickness (x), and 2nm gate thickness (+). the poor subthreshold slope, and may be reduced by increasing the length of the intrinsic region from the value of 7nm used here. In previous works, both experimental [64] and theoretical [65], the intrinsic regions were at least 20 nm in length. This work indicates that there is a lower limit to the channel length in order to realize the benefit of an exceptional subthreshold slope, i.e., using a p-i-n device does not necessarily imply that this figure-of-merit will be better than the thermionic limit. From the figure, it is also clear that neglecting the azimuthal potential variation leads to a fairly small error in the D G terminal characteristic. In addition, note that the gate thickness has a similarly small effect for the thicknesses considered here. For the SC case, however, the gate thickness plays a larger role, with roughly 1.6 times as much current if a zero gate thickness is assumed. This discrepancy exists even for the aggressive value of 2 nm for the gate thickness simulated here, which would lead to excessive series resistance in a practical device. In order to understand this increase, consider the energy band diagrams for the semi-cylindrical geometries in Fig. 5.12. The solid line is for a gate thickness of 2 nm, and it may be seen that the inter-band tunneling distance is greater 83 -10 0 10 20 Distance from source (nm) Figure 5.12: Energy band diagrams for the semi-cylindrical-gate geometry with a gate thickness of 2nm (solid) and Onm (dashed). in this case than for the zero-gate-thickness device, shown by the dashed line. This is due to the fringing field from the gate, which extends to greater distances in the thick-gate device. The steeper band-bending in the thin-gate case results in increased tunneling, and hence, more current. Turning now to the atomic current, two of the DG devices are considered. For the case where the axial potential is used, the current shows no azimuthal variation, and is distributed evenly around the circumference of the CN. For the device with zero gate thickness, where the proper surface potential at each atomic site was utilized, considerable structure emerges, as shown in Fig. 5.13. While the ultimate utility of this calculation is not yet known, perhaps this approach could be used to investigate the effect of defects on CNFET behaviour, i.e., the impact of defects on the electron flow pattern through the CN, or the flow patterns across junctions between CNs of differing chirality. Finally, the NEGF model is used to compare to results generated3 from using the Flietner dispersion relation throughout the full energy range [76] on a p-i-n device. A cylindrical structure is used for this comparison, such that the ID Schrodinger solution used in Ref. [76] applies. The device consists of a 7 nm intrinsic channel connected to doped CN leads, and a 2 nm Si02 layer between 3Doped-contact results using the Flietner dispersion relation were provided by L . Castro of the University of Bri t ish Columbia. 84 1 x 10 Mrad) 0 0 z ( m ) Figure 5.13: Inter-atomic current for a double-gate device with zero gate thickness. the channel and the gate. The doping density is 0.1 nm - 2 as for the previous NEGF results. As Fig. 5.14 shows, the Flietner relation provides very good agreement with the NEGF results in the low-current regime. At higher current, the discrepancy between the predictions grows to a maximum value of about 32% at VGS — —0.06 V, where it is noted that the Flietner band structure results in an overestimate of the current. It remains to be seen how this comparison will turn out when considering non-cylindrical geometries. 85 10 ~ i o " 6 10 - 0 . 2 o NEGF ° Flietner 0.2 0.4 0.6 0 .8 Figure 5.14: Comparison of NEGF results (o) to those using the Flietner dispersion relation (•) for a cylindrical p-i-n device. 86 Chapter 6 Conclusions In this thesis, we have presented a variety of CNFET models valid in different operating regimes, and subject to successively less restrictive approximations. We have focused on including quantum phenomena that had been neglected at the time this work commenced. For example, the "top-of-the-barrier" [23,52] model, while it may be useful for negative barrier devices, it does not deal with tunneling or quantum resonance. Also, the quasi-Fermi level model [68,69,71], while it was one of the earliest models to predict bipolar conduction at high Vos, predicted a much earlier onset of this behaviour than the more rigorous effective-mass model due to an unphysical thinning of the Schottky barriers at high bias. 6.1 Current Progress Our first contribution was to provide equilibrium results, using an analytical expression for the Poisson solution, from which Castro et ai, derived a compact model [51] that incorporated the WKB approximation to tunneling. In Ref. [77], we used this model to present some preliminary current-voltage data, and also showed the dependence of the equilibrium band structure on oxide permittivity and contact workfunction. A modified WKB result was derived to also account for above-barrier reflection, something that a classical treatment would also ignore, and this led to significant improvements in Castro's compact model [53]. Following this, we proceeded to develop a coupled solution to the Poisson and effective-mass Schrodinger equations [72]. This work served as the early benchmark against which compact models could be validated, and also served as the method for performing small-signal analyses [73-75,103] and switching-time calculations [55]. In addition, we derived, in detail, expressions for the so-called quantum capacitance [89], which had been previously used in qualitative discussions of CNFETs [23,30,90]. In so doing, we observed that this capacitance is not truly quantized due to the form 87 of the one-dimensional DOS. Our expression may be used to assist in the quantitative modeling of CNFETs. In order to deal with the non-parabolicity of the CN bands, researchers started to investigate simulations using the NEGF formalism. In Ref. [56], a method was presented for simulating devices where the potential on the surface of the CN was not a function of azimuthal angle. In our work, we extended that technique by considering a real-space basis where the potential at each atom was computed via a three-dimensional finite element solution. While completing this method, we became aware of Refs. [57,58] where essentially the same technique has been used. One difference in our implementations was in the computation of the self-energy matrices, where we derived a quadratic matrix equation [117] which may be solved with a generalization of Newton's method. In Refs. [56-58], a recursive algorithm was employed [118], which we have found to be less computationally efficient than our method. In addition, these works have used a Gummel technique in order to obtain convergence, whereas we have chosen an approximate Jacobian [132] in Newton's method, which seems to be more stable, albeit more time-consuming, in the cases that we have investigated thus far. Despite the similarities in algorithm, we have further expanded on the work of Fiori et al, where devices were simulated assuming that the contact thickness played no role in the device characteristics. This was implicitly assumed by the imposition of Neumann boundary conditions along surfaces that extended from the nearest face of the gate-contact to the CN channel [57,58]. In this work, we investigated the effect of contact thickness, and also used our finite element Poisson solver in order to consider curved gate contacts. In addition, Fiori et al., present results only for the (11,0) CN, with a diameter of 0.86 nm. It is known that CNs with diameters below lnm undergo curvature-induced hybridization of the a and IT states [137] which significantly alter the CN band structure from the result of considering only the TT states. In most of the literature, and in this thesis, we consider only the n states. As a result, we present results primarily for the (16,0) CN, with a diameter of 1.25 nm, for which our calculations should be reasonable. Finally, while the ultimate utility of this aspect of the thesis is not yet clear, we have also derived the current operator for a carbon nanotube such that we may show the current between individual atoms of the CN. From these contributions to the understanding of CNFET behaviour, we may conclude the following: 88 1. equilibrium calculations may be performed to gain insight into CNFET operation, and to aid in the development of compact models; 2. the tight-binding approach may be used to various levels of approximation in order to simulate CNFETs to varying degrees of physical rigour and computational effort; 3. effective-mass Schrodinger-Poisson solutions may be used for DC, small-signal AC, and switch-ing time calculations; 4. an integral expression may be derived for the quantum capacitance; 5. significant optimization of CNFET geometries will be required in order to obtain switching times that rival existing Si CMOS, such as, for example, using parallel CNs in order to increase the drive current, or decreasing the source and drain contact thicknesses in order to reduce the extrinsic capacitances; 6. for large-signal behaviour, it is desirable to have a large bandgap, in order to obtain a good 2 O N / I O F P , and a thin oxide, in order to decrease the switching time by increasing X O N ; 7. the non-equilibrium Green's function formalism may be used with an atomistic, real-space, Hamiltonian in order to treat multi-band effects, and azimuthal asymmetries in CNFETs; 8. a quadratic matrix equation may be solved in order to improve the computation time of the self-energy matrices representing the semi-infinite leads connected to the CN channel; 9. contact thickness may play an important role in a device's behaviour for geometries where the gate-channel coupling is strong; 10. the current operator for a CN may be used to investigate electron flow patterns at the inter-atomic level; 11. while the azimuthal potential variation does not significantly affect the terminal characteristics, for the devices considered here, the inter-atomic current shows considerable structure, which would be missed if azimuthal effects are ignored; 12. in order to realize the exceptional subthreshold slope made possible by the use of a p-i-n structure, the channel length needs to be long enough to suppress the significant inter-band tunneling of carriers in the OFF-state; 89 13. use of the Flietner dispersion relation for a cylindrical device results in reasonable agreement with the NEGF results, particularly in the low-current regime, with greater error when the device is ON. 6.2 Future Work Using these advances as the foundation for further study, we finally consider some possibilities for future work. Of course, there is initial work on geometry optimization that may be performed directly from this work in order to confirm, and extend upon, the results of Refs. [109,138]. Along with this, a detailed evaluation of existing models may be performed in order to validate less computationally intensive methods than the thorough NEGF algorithm. This work was begun in Chapter 5, where we compared to results generated using the Flietner dispersion relation, and may be continued to include a larger bias range, and non-cylindrical devices. Most importantly, these models need to be rigorously validated against experimental results. This has been impeded by the many technological challenges involved with CNFET fabrication; however, our atomistic approach will better lend itself to the simulation of experimental structures than the earlier methods which relied on cylindrical symmetry. Models having the same general form as employed here have been compared, and fit, to individual device characteristics [56,109,112,139], and this is an encouraging sign of the validity and usefulness of our modeling approach. It also gives confidence in the quantitatively predictive quality of our model, but the ultimate validation of this will require a systematic experimental investigation which has not yet been achieved. In terms of continuing the progression of these models towards more physically rigorous sit-uations, the inclusion of electron-phonon interactions has started to be considered by some au-thors [12,65], and this may also be incorporated into our model through additional self-energy matrices [12,94]. In the DC case, this is not expected to impact on performance, but will have an effect on the AC behaviour [12].' Another extension, already alluded to in this work, would be to consider the role of defects in CNs, and to treat junctions between dissimilar tubes. Such defects have been shown to produce single-electron transistor characteristics [140] due to local band structure modifications [141], and have also been predicted to induce metal-semiconductor transitions in some tubes [141,142]. This may assist with the integration of devices if a single tube could be used to realize a metal-contacted 90 ? semiconductor channel [142]. Knowledge of the carrier flow patterns may be helpful in this investi-gation. Deformations of pristine nanotubes due to their physical interaction with a substrate could also be treated by modifying the overlap integral 7 such that it varies with azimuthal angle. A density functional theory calculation could be used in order to obtain appropriate values. Further, while DC simulations, such as those described herein, can be used to perform small-signal and switching-time analyses, it would be helpful to perform time-dependent simulations in order to validate previous results [73-76,103]. For a switching calculation, for example, our DC solution could be used as an initial condition, with the evolution of the wavefunction given by the time-dependent Schrodinger equation. Finally, it may be important to consider electron-electron interactions, as may lead to Liittinger liquid behaviour [143]; excitonic interactions, as may be important in electro-optical devices; or even just to include the effects of exchange and correlation, as in a density functional theory approach. Some investigations into the Liittinger liquid phenomena have been performed, e.g., Refs. [144,145], but not much in the context of CNFETs except perhaps for a qualitative treatment in Ref. [146]. The importance of this phenomena for CNFETs is not yet clear as it does not appear to be present at some energy scales [147], and electron-electron interactions are expected to have little effect on conductance when the contacts are highly transparent [148], such as would be desirable in a CNFET. However, many-body effects have also been shown to significantly modify the CN bandgap [149], and if this occurs in CNFETs, it would be critical to incorporate this effect. The solid framework provided by this thesis will be an excellent starting point for these, and possibly other, further investigations. By continuing to successively strip away some of the more re-strictive assumptions of our model, significant progress may continue towards a better understanding of these fascinating devices. 91 Bibliography [1] Semiconductor Industry Association. International Technology Roadmap for Semiconductors—2004 Update: Overview and Summaries, 2004. [Online.] Available: http://www.itrs.net/Common/2004Update/2004Update.htm. [2] A. Oberlin, M. Endo, and A. T. Koyama. Filamentous growth of carbon through benzene decomposition. J. Cryst. Grow., 32(3):335-349, 1976. [3] Sumio Iijima. Helical microtubules of graphitic carbon. Nature, 354:56-58, 1991. [4] R. Saito, G. Dresselhaus, and M. S.1 Dresselhaus. Physical Properties of Carbon Nanotubes. Imperial College Press, London, 1st edition, 1998. [5] M. S. Dresselhaus, G. Dresselhaus, and P. C. Eklund. Science of Fullerenes and Carbon Nanotubes. Academic Press, Toronto, 1996. [6] Keivan Esfarjani, Amir A. Farajian, Yuichi Hashi, and Yoshiyuki Kawazoe. Electronic, trans-port and mechanical properties of carbon nanotubes. In Y. Kawazoe, T. Kondow, and K. Ohno, editors, Clusters and Nanomaterials—Theory and Experiment, Springer Series in Cluster Physics, pages 187-220. Springer-Verlag, Berlin, 2002. [7] Sander J. Tans, Alwin R. M. Verschueren, and Cees Dekker. Room-temperature transistor based on a single carbon nanotube. Nature, 393:49-52, 1998. [8] Paul L. McEuen, Michael S. Fuhrer, and Hongkun Park. Single-walled carbon nanotube electronics. IEEE Trans. Nanotechnol, l(l):78-85, 2002. [9] T. Durkop, E. Cobas, and M. S. Fuhrer. High-mobility semiconducting nanotubes. In Molecular Nanostructures, pages 524-527, 2003. 92 [10] Ali Javey, Jing Guo, Magnus Paulsson, Qian Wang, David Mann, Mark Lundstrom, and Hongjie Dai. High-field, quasi-ballistic transport in short carbon nanotubes. Phys. Rev. Lett., 92(10):106804-l-106804-4, 2004. [11] Yu-Ming Lin, Joerg Appenzeller, Zhihong Chen, Zhi-Gang Chen, Hui-Ming Cheng, and Phae-don Avouris. High-performance dual-gate carbon nanotube FETs with 40-nm gate length. IEEE Electron Device Lett, 26(ll):823-825, 2005. [12] Jing Guo. A quantum-mechanical treatment of phonon scattering in carbon nanotube tran-sistors. J. Appl. Phys., 98:063519-1-063519-6, 2005. [13] Akin Akturk, Gary Pennington, and Neil Goldsman. Quantum modeling and proposed designs of CNT-embedded nanoscale MOSFETs. IEEE Trans. Electron Devices, 52(4):577-584, 2005. [14] Khairul Alam and Roger K. Lake. Leakage and performance of zero-Schottky-barrier carbon nanotube transistors. J. Appl. Phys., 98:064307-1-064307-8, 2005. [15] Vasili Perebeinos, J. Tersoff, and Phaedon Avouris. Mobility in semiconducting carbon nan-otubes at finite carrier density. Nano Lett., 6(2):205-208, 2006. [16] David L. Pulfrey and N. Garry Tarr. Introduction to Microelectronic Devices. Prentice Hall Series in Solid State Physical Electronics. Prentice Hall, New Jersey, 1989. [17] Yi Cui, Zhaohui Zhong, Deli Wang, Wayne U. Wang, and Charles M. Lieber. High performance silicon nanowire field effect transistors. Nano Lett., 3(2):149-152, 2003. [18] Cees Dekker. Carbon nanotubes as molecular quantum wires. Phys. Today, 52(5):22-28, 1999. [19] Zhen Yao, Charles L. Kane, and Cees Dekker. High-field electrical transport in single-wall carbon nanotubes. Phys. Rev. Lett., 84(13):2941-2944, 2000. [20] A. Tilke, L. Pescini, A. Erbe, H. Lorenz, and R. H. Blick. Electron-phonon interaction in suspended highly doped silicon nanowires. Nanotechnology, 13:491-494, 2002. [21] Sang-Hyun Oh, Don Monroe, and J. M. Hergenrother. Analytic description of short-channel effects in fully-depleted double-gate and cylindrical, surrounding-gate MOSFETs. IEEE Elec-tron Device Lett, 21(9):445-447, 2000. 93 [22] Brian Winstead and Umberto Ravaioli. Simulation of Schottky barrier MOSFET's with a coupled quantum injection/Monte Carlo technique. IEEE Trans. Electron Devices, 47(6): 1241-1246, 2000. [23] Jing Guo, Sebastien Goasguen, Mark Lundstrom, and Supriyo Datta. Metal-insulator-semiconductor electrostatics of carbon nanotubes. Appl. Phys. Lett., 81(8):1486-1488, 2002. [24] Wolfgang Hoenlein, Franz Kreupl, Georg Stefan Duesberg, Andrew Peter Graham, Maik Liebau, Robert Viktor Seidel, and Eugen Unger. Carbon nanotube applications in micro-electronics. IEEE Trans. Compon. Pack. T., 27(4):629-634, 2004. [25] Hou T. Ng, J. Han, Toshishige Yamada, P. Nguyen, Yi P. Chen, and M. Meyyappan. Single crystal nanowire vertical surround-gate field-effect transistor. Nano Lett, 4(7): 1247-1252, 2004. [26] Ali Javey, Hyoungsub Kim, Markus Brink, Qian Wang, Ant Ural, Jing Guo, Paul Mclntyre, Paul McEuen, Mark Lundstrom, and Hongjie Dai. High-/c dielectrics for advanced carbon-nanotube transistors and logic gates. Nature Mater., 1:241-246, 2002. [27] B. M. Kim, T. Brintlinger, E. Cobas, M. S. Fuhrer, Haimei Zheng, Z. Yu, R. Droopad, J. Ram-dani, and K. Eisenbeiser. High-performance carbon nanotube transistors on SrTiC-3/Si sub-strates. Appl. Phys. Lett, 84(11):1946-1948, 2004. [28] Sami Rosenblatt, Yuval Yaish, Jiwoong Park, Jeff Gore, Vera Sazonova, and Paul L. McEuen. High performance electrolyte gated carbon nanotube transistors. Nano Lett, 2(8):869-872, 2002. [29] S. B. Cronin, R. Barnett, M. Tinkham, S. G. Chou, O. Rabin, M. S. Dresselhaus, A. K. Swan, M. S. Unlu, and B. B. Goldberg. Electrochemical gating of individual single-wall carbon nanotubes observed by electron transport measurements and resonant Raman spectroscopy. Appl. Phys. Lett, 84(12):2052-2054, 2004. [30] Anisur Rahman, Jing Guo, Supriyo Datta, and Mark S. Lundstrom. Theory of ballistic nan-otransistors. IEEE Trans. Electron Devices, 50(9):1853-1864, 2003. 94 [31] Ming Zheng, Anand Jagota, Ellen D. Semke, Bruce A. Diner, Robert S. McLean, Steve R. Lustig, Raymond E. Richardson, and Nancy G. Tassi. DNA-assisted dispersion and separation of carbon nanotubes. Nature Mater., 2:338-342, 2003. [32] Gang Lu, Paul Maragakis, and Efthimios Kaxiras. Carbon nanotube interaction with DNA. Nano Lett, 5(5):897-900, 2005. [33] Jing Kong, Nathan R. Franklin, Chongwu Zhou, Michael G. Chapline, Shu Peng, Kyeongjae Cho, and Hongjie Dai. Nanotube molecular wires as chemical sensors. Science, 287:622-625, 2000. [34] Keith Bradley, Jean-Christdphe P. Gabriel, Alexander Star, and George Griiner. Short-channel effects in contact-passivated nanotube chemical sensors. Appl. Phys. Lett., 83(18):3821-3823, 2003. [35] Philip G. Collins and Phaedon Avouris. Nanotubes for electronics. Sci. Am., pages 62-69, December 2000. [36] R. Gutierrez, G. Fagas, G. Cuniberti, F. Grossmann, R. Schmidt, and K. Richter. Theory of an all-carbon molecular switch. Phys. Rev. B, 65:113410-1-113410-4, 2002. [37] R. Martel, T. Schmidt, H. R. Shea, T. Hertel, and Ph. Avouris. Single- and multi-wall carbon nanotube field-effect transistors. Appl. Phys. Lett, 73(17):2447-2449, 1998. [38] Hongjie Dai. Nanotube growth and characterization. In Mildred S. Dresselhaus, Gene Dres-selhaus, and Phaedon Avouris, editors, Carbon Nanotubes: Synthesis, Structure, Properties, and Applications, volume 80 of Topics in Applied Physics, pages 29-53. Springer, Berlin, 2001. [39] Hongjie Dai. Carbon nanotubes: Synthesis, integration, and properties. Acc. Chem. Res., 35(12):1035-1044, 2002. [40] Ralph Krupke, Frank Hennrich, Hilbert v. Lohneysen, and Manfred M. Kappes. Separation of metallic from semiconducting single-walled carbon nanotubes. Science, 301:344-347, 2003. [41] Michael S. Strano, Christopher A. Dyke, Monica L. Usrey, Paul W. Barone, Matthew J. Allen, Hongwei Shan, Carter Kittrell, Robert H. Hauge, James M. Tour, and Richard E. Smal-95 ley. Electronic structure control of single-walled carbon nanotube functionalization. Science, 301:1519-1522, 2003. [42] Yutaka Ohno, Shigeru Kishimoto, Takashi Mizutani, Toshiya Okazaki, and Hisanori Shinohara. Chirality assignment of individual single-walled carbon nanotubes in carbon nanotube field-effect transistors by micro-photocurrent spectroscopy. Appl. Phys. Lett., 84(8): 1368-1370, 2004. [43] Aleksey N. Kolmogorov, Vincent H. Crespi, Monica H. Schleier-Smith, and James C. Ellenbo-gen. Nanotube-substrate interactions: Distinguishing carbon nanotubes by the helical angle. Phys. Rev. Lett., 92(8):085503-l-085503-4, 2004. [44] Yiming Li, David Mann, Marco Rolandi, Woong Kim, Ant Ural, Steven Hung, Ali Javey, Jien Cao, Dunwei Wang, Erhan Yenilmez, Qian Wang, James F. Gibbons, Yoshio Nishi, and Hongjie Dai. Preferential growth of semiconducting single-walled carbon nanotubes by a plasma enhanced CVD method. Nano Lett, 4(2):317-321, 2004. [45] J. Appenzeller, J. Knoch, R. Martel, V. Derycke, S. Wind, and Ph. Avouris. Short-channel like effects in Schottky barrier carbon nanotube field-effect transistors. In IEDM Tech. Digest, pages 285-288, 2002. [46] Ph. Avouris, R. Martel, S. Heinze, M. Radosavljevic, S. Wind, V. Derycke, J. Appenzeller, and J. Terso. The role of Schottky barriers on the behavior of carbon nanotube field-effect transistors. In H. Kuzmany, J. Fink, M. Mehring, and S. Roth, editors, Structural and Elec-tronic Properties of Molecular. Nano structures, XVI International Winterschool on Electronic Properties of Novel Materials, pages 508-512, 2002. [47] S. Heinze, J. Tersoff, R. Martel, V. Derycke, J. Appenzeller, and Ph. Avouris. Carbon nan-otubes as Schottky barrier transistors. Phys. Rev. Lett, 89(10):106801-l-106801-4, 2002. [48] J. Knoch and J. Appenzeller. Impact of the channel thickness on the performance of Schottky barrier metal-oxide-semiconductor field-effect transistors. Appl. Phys. Lett., 81(16):3082-3084, 2002. 96 [49] J. Appenzeller, J. Knoch, and Ph. Avouris. Carbon nanotube field-effect transistors—an example of an ultra-thin body Schottky barrier device. In Proc. TMS 61st Annual Device Research Conference, pages 167-170, Salt Lake City, U.S.A., June 2003. [50] S. Heinze, M. Radosavljevic, J. Tersoff, and Ph. Avouris. Unexpected scaling of the perfor-mance of carbon nanotube Schottky-barrier transistors. Phys. Rev. B, 68:235418-1-235418-5, 2003. [51] L. C. Castro, D. L. John, and D. L. Pulfrey. Towards a compact model for Schottky-barrier nanotube FETs. In Proc. IEEE Conf. on Optoelectronic and Microelectronic Materials and Devices, pages 303-306, Sydney, Australia, December 2002. [52] Jing Guo, Mark Lundstrom, and Supriyo Datta. Performance projections for ballistic carbon nanotube field-effect transistors. Appl. Phys. Lett, 80(17):3192-3194, 2002. [53] L. C. Castro, D. L. John, and D. L. Pulfrey. An improved evaluation of the DC performance of carbon nanotube field-effect transistors. Smart Mater. Struct, 15(1):S9-S13, 2006. [54] M. P. Anantram, S. Datta, and Yongqiang Xue. Coupling of carbon nanotubes to metallic contacts. Phys. Rev. B, 61(20):14219-14225, 2000. [55] D. L. John and D. L. Pulfrey. Switching-speed calculations for Schottky-barrier carbon nan-otube field-effect transistors. J. Vac. Sci. TechnoL A, 24(3):708-712„ 2006. [56] Jing Guo. Carbon nanotube electronics: Modeling, physics, and applications. PhD thesis, Purdue University, West Lafayette, U.S.A., 2004. [57] G. Fiori, G. Iannaccone, and G. Klimeck. Performance of carbon nanotube field effect transis-tors with doped source and drain extensions and arbitrary geometry. In IEDM Tech. Digest, pages 529-532, December 2005. [58] Gianluca Fiori, Giuseppe Iannaccone, Mark Lundstrom, and Gerhard Klimeck. Three-dimensional atomistic simulation of carbon nanotube FETs with realistic geometry. In Proc. European Solid-State Device Research Conference, Grenoble, France, September 2005. [59] M. Radosavljevic, J. Appenzeller, Ph. Avouris, and J. Knoch. High performance of potassium n-doped carbon nanotube field-effect transistors. Appl. Phys. Lett., 84(18):3693-3695, 2004. 97 [60] Jia Chen, Christian Klinke, Ali Afzali, Kevin Chan, and Phaedon Avouris. Self-aligned carbon nanotube transistors with novel chemical doping. In IEDM Tech. Digest, pages 695-698, December 2004. [61] Yu-Ming Lin, Joerg Appenzeller, and Phaedon Avouris. Novel carbon nanotube FET design with tunable polarity. In IEDM Tech. Digest, pages 687-690, December 2004. [62] Ali Javey, Ryan Tu, Damon B. Farmer, Jing Guo, Roy G. Gordon, and Hongjie Dai. High performance n-type carbon nanotube field-effect transistors with chemically doped contacts. Nano Lett, 5(2):345-348, 2005. [63] J. Appenzeller, Y.-M. Lin, J. Knoch, and Ph. Avouris. Band-to-band tunneling in carbon nanotube field-effect transistors. Phys. Rev. Lett., 93(19):196805-l-196805-4, 2004. [64] Joerg Appenzeller, Yu-Ming Lin, Joachim Knoch, Zhihong Chen, and Phaedon Avouris. Com-paring carbon nanotube transistors-the ideal choice: A novel tunneling device design. IEEE Trans. Electron Devices, 52(12):2568-2576, 2005. [65] Siyuranga O. Koswatta, Dmitri E. Nikonov, and Mark S. Lundstrom. Computational study of carbon nanotube p-i-n tunnel FETs. In IEDM Tech. Digest, pages 525-528, December 2005. [66] Neil W. Ashcroft and N. David Mermin. Solid State Physics. Harcourt College Publishers, New York, 1st edition, 1976. [67] Jeroen W. G. Wildoer, Liesbeth C. Venema, Andrew G. Rinzler, Richard E. Smalley, and Cees Dekker. Electronic structure of atomically resolved carbon nanotubes. Nature, 391:59-62, 1998. [68] Jason Clifford, D. L. John, and David L. Pulfrey. Bipolar conduction and drain-induced barrier thinning in carbon nanotube FETs. IEEE Trans. Nanotechnol, 2(3):181-185, 2003. [69] Jason Paul Clifford. A self-consistent numerical model for bipolar transport in carbon nanotube field-effect transistors. Master's thesis, University of British Columbia, Vancouver, Canada, December 2003. [70] L. C. Castro, D. L. John, and D. L. Pulfrey. Carbon nanotube transistors: An evaluation. In Proc. SPIE Conf. Device and Process Technologies for MEMS, Microelectronics, and Photonics III, volume 5276, pages 1-10, Perth, Australia, April 2004. 98 [71] J. P. Clifford, D. L. John, L. C. Castro, and D. L. Pulfrey. Electrostatics of partially gated carbon nanotube FETs. IEEE Trans. Nanotechnol, 3(2):281-286, 2004. [72] D. L. John, L. C. Castro, P. J. S. Pereira, and D. L. Pulfrey. A Schrodinger-Poisson solver for modeling carbon nanotube FETs. In Tech. Proc. of the 2004 NSTI Nanotechnology Conf. and Trade Show, volume 3, pages 65-68, Boston, U.S.A., March 2004. [73] L. C. Castro, D. L. John, D. L. Pulfrey, M. Pourfath, A. Gehring, and H. Kosina. Method for predicting fT for carbon nanotube FETs. IEEE Trans. Nanotechnol., 4(6):699-704, 2005. [74] L. C. Castro and D. L. Pulfrey. Extrapolated / m a x for carbon nanotube field-effect transistors. Nanotechnology, 17(l):300-304, 2006. [75] L. C. Castro, D. L. Pulfrey, and D. L. John. High-frequency capability of Schottky-barrier carbon nanotube FETs. Solid-State Phenomena, 2005. Accepted December 22, 2005. [76] Leonardo C. Castro. Modeling of carbon nanotube field-effect transistors. PhD thesis, Univer-sity of British Columbia, Vancouver, Canada, 2006. [77] D. L. John, Leonardo C. Castro, Jason Clifford, and David L. Pulfrey. Electrostatics of coaxial Schottky-barrier nanotube field-effect transistors. IEEE Trans. Nanotechnol, 2(3):175-180, 2003. [78] Christopher P. Auth and James D. Plummer. Scaling theory for cylindrical, fully-depleted, surrounding-gate MOSFET's. IEEE Electron Device Lett, 18(2):74-76, 1997. [79] Jing Guo, Ali Javey, Hongjie Dai, Supriyo Datta, and Mark Lundstrom. Predicted performance advantages of carbon nanotube transistors with doped nanotubes as source/drain. [Online.] Available: http://arxiv.org/pdf/cond-mat/0309039, 2003. [80] Francois Leonard and J. Tersoff. Role of Fermi-level pinning in nanotube Schottky diodes. Phys. Rev. Lett, 84(20):4693-4696, 2000. [81] Francois Leonard and J. Tersoff. Dielectric response of semiconducting carbon nanotubes. Appl. Phys. Lett., 81(25):4835-4837, 2002. 99 [82] John David Jackson. Classical Electrodynamics. John Wiley and Sons, Toronto, 3rd edition, 1999. [83] D. L. McGuire and D. L. Pulfrey. Error analysis of boundary condition approximations in the modeling of coaxially-gated carbon nanotube field-effect transistors. Phys. Stat. Sol. (a), 2006. Accepted February 21, 2006. [84] E. M. Freeman and D. A. Lowther. A novel mapping technique for open boundary finite element solutions to Poisson's equation. IEEE Trans. Magn., 24(6):2934-2936, 1988. [85] FEMLAB reference manual, http://www.femlab.com, 2002. [86] Ali Javey, Moonsub Shim, and Hongjie Dai. Electrical properties and devices of large-diameter single-walled carbon nanotubes. Appl. Phys. Lett, 80(6): 1064-1066, 2002. [87] Takeshi Nakanishi, Adrian Ba,chtold, and Cees Dekker. Transport through the interface be-tween a semiconducting carbon nanotube and a metal electrode. Phys. Rev. B, 66:073307-1-073307-4, 2002. [88] Serge Luryi. Quantum capacitance devices. Appl. Phys. Lett., 52(6):501-503, 1988. [89] D. L. John, L. C. Castro, and D. L. Pulfrey. Quantum capacitance in nanoscale device mod-eling. J. Appl. Phys., 96(9):5180-5184, 2004. [90] P. J. Burke. An RF circuit model for carbon nanotubes. IEEE Trans. Nanotechnol., 2(l):55-58, 2003. [91] Peter J. Burke. AC performance of nanoelectronics: Towards a ballistic THz nanotube tran-sistor. Solid-State Electron., 48:1981-1986, 2004. [92] D. L. Pulfrey. Modeling high-performance HBTs. In Patrick Roblin and Hans Rohdin, editors, High-Speed Heterostructure Devices, chapter 18. Cambridge University Press, 2002. [93] David K. Ferry and Stephen M. Goodnick. Transport in Nanostructures. Cambridge University Press, New York, 1997. 100 [94] Supriyo Datta. Electronic Transport in Mesoscopic Systems, volume 3 of Cambridge Studies in Semiconductor Physics and Microelectronic Engineering. Cambridge University Press, New York, 1995. [95] M. H. Holmes. Introduction to Perturbation Methods, volume 20 of Texts in Applied Mathe-matics. Springer-Verlag, New York, 1995. [96] David J. Griffiths. Introduction to Quantum Mechanics. Prentice Hall, New Jersey, 1994. [97] David Yuk Kei Ko and J. C. Inkson. Matrix method for tunneling in heterostructures: Reso-nant tunneling in multilayer systems. Phys. Rev. B, 38(14):9945-9951, 1988. [98] Wayne W. Lui and Masao Fukuma. Exact solution of the Schrodinger equation across an arbitrary one-dimensional piecewise-linear potential barrier. J. Appl. Phys., 60(5): 1555-1559, 1986. [99] J. Tersoff. Schottky barrier heights and the continuum of gap states. Phys. Rev. Lett., 52(6):465-468, 1984. [100] Lee W. Johnson and R. Dean Riess. Numerical Analysis. Addision-Wesley, Don Mills, Ontario, 1977. [101] Zhen Yao, Cees Dekker, and Phaedon Avouris. Electrical transport through single-wall car-bon nanotubes. In Mildred S. Dresselhaus, Gene Dresselhaus, and Phaedon Avouris, editors, Carbon Nanotubes, volume 80 of Topics Appl. Phys., pages 147-171. Springer-Verlag, Berlin, 2001. [102] Hang Si. TetGen: A quality tetrahedral mesh generator and three-dimensional Delaunay tri-angulator (Version 1.3), June 2005. [Online.] Available: http://tetgen.berlios.de. [103] D. L. Pulfrey, D. L. John, and L. C. Castro. Carbon nanotube field-effect transistors: AC performance capabilities. In Proc. 13th Intl. Workshop Phys. Semiconductor Dev., pages 7-13, Delhi, India, December 2005. Invited. [104] Richard Martel, Hon-Sum Philip Wong, Kevin Chan, and Phaedon Avouris. Carbon nanotube field effect transistors for logic applications. In IEDM Tech. Digest, pages 159-162, 2001. 101 [105] R. Martel, V. Derycke, C. Lavoie, J. Appenzeller, K. K. Chan, J. Tersoff, and Ph. Avouris. Ambipolar electrical transport in semiconducting single-wall carbon nanotubes. Phys. Rev. Lett, 87(25):256805-l-256805-4, 2001. [106] R. Martel, V. Derycke, J. Appenzeller, S. Wind, and Ph. Avouris. Carbon nanotube field-effect transistors and logic circuits. In Proc. Design Automation Conf., pages 94-98, June 2002. [107] M. Radosavljevic, S. Heinze, J. Tersoff, and Ph. Avouris. Thermal carrier injection into ambipolar carbon nanotube field effect transistors. In Molecular Nanostructures, March 2003. [108] Robert Chau, Suman Datta, Mark Doczy, Brian Doyle, Ben Jin, Jack Kavalieros, Amlan Majumdar, Matthew Metz, and Marko Radosavljevic. Benchmarking nanotechnology for high-performance and low-power logic transistor applications. IEEE Trans. Nanotechnol, 4(2): 153-158, 2005. [109] Jing Guo, Ali Javey, Hongjie Dai, and Mark Lundstrom. Performance analysis and design optimization of near ballistic carbon nanotube field-effect transistors. In IEDM Tech. Digest, pages 703-706, December 2004. [110] Ali Javey, Jing Guo, Damon B. Farmer, Qian Wang, Erhan Yenilmez, Roy G. Gordon, Mark Lundstrom, and Hongjie Dai. Self-aligned ballistic molecular transistors and electrically par-allel nanotube arrays. Nano Lett, 4(7):1319-1322, 2004. [Ill] Jijun Zhao, Jie Han, and Jian Ping Lu. Work functions of pristine and alkali-metal intercalated carbon nanotubes and bundles. Phys. Rev. B, 65:193401-1-193401-4, 2002. [112] Zhihong Chen, Joerg Appenzeller, Joachim Knoch, Yu-Ming Lin, and Phaedon Avouris. The role of metal-nanotube contact in the performance of carbon nanotube field-effect transistors. Nano Lett, 5(7):1497-1502, 2005. [113] Jing Guo. Personal communication. [114] Supriyo Datta. Nanoscale device modeling: The Green's function method. Superlattices Mi-crostruct, 28(4):253-278, 2000. [115] H. Flietner. The E(k) relation for a two-band scheme of semiconductors and the application to the metal-semiconductor contact. Phys. Stat Sol. (b), 54:201-208, 1972. [116] A. Persson and R. M. Cohen. Reformulated Hamiltonian for nonparabolic bands in semicon-ductor quantum wells. Phys. Rev. B, 38(8):5568-5575, 1988. [117] D. L. John and D. L. Pulfrey. Green's function calculations for semi-infinite carbon nanotubes. Phys. Stat. Sol. (b), 243(2):442-448, 2006. [118] A. Svizhenko, M. P. Anantram, T. R. Govindan, B. Biegel, and R. Venugopal. Two-dimensional quantum mechanical modeling of nanotransistors. J. Appl. Phys., 91(4):2343-2354, 2002. [119] E. N. Economou. Green's Functions in Quantum Physics, volume 7 of Springer Series in Solid-State Sciences. Springer-Verlag, New York, 2nd edition, 1983. [120] Erich Zauderer. Partial Differential Equations of Applied Mathematics. John Wiley and Sons, Toronto, 2nd edition, 1998. [121] E. B. Saff and A. D. Snider. Fundamentals of Complex Analysis for Mathematics, Science, and Engineering. Prentice Hall, New Jersey, 2nd edition, 1993. [122] Ira N. Levine. Quantum Chemistry. Prentice Hall, New Jersey, 5th edition, 2000. [123] Ian Appelbaum, Tairan Wang, J. D. Joannopoulos, and V. Narayanamurti. Ballistic hot-electron transport in nanoscale semiconductor heterostructures: Exact self-energy of a three-dimensional periodic tight-binding Hamiltonian. Phys. Rev. B, 69:165301-1-165301-6, 2004. [124] R. Venugopal, Z. Ren, S. Datta, M. S. Lundstrom, and D. Jovanovic. Simulating quan-tum transport in nanoscale transistors: Real versus mode-space approaches. J. Appl. Phys., 92(7):3730-3739, 2002. [125] Nicholas J. Higham and Hyun-Min Kim. Solving a quadratic matrix equation by Newton's method with exact line searches. SI AM J. Matrix Anal. Appl, 23(2):303-316, 2001. [126] Judith D. Gardiner, Alan J. Laub, James J. Amato, and Cleve B. Moler. Solution of the Sylvester matrix equation axb1 + cxd* = e. ACM Trans. Math. Softw., 18(2):223-231, 1992. [127] J. Knoch and J. Appenzeller. A novel concept for field-effect transistors—the tunneling carbon nanotube FET. In Proc. TMS 63rd Annual Device Research Conference, pages 100-103, Santa Barbara, U.S.A., June 2005. 103 [128] Albert Messiah. Quantum Mechanics. Dover, New York, 1999. [129] Frederick W. Byron Jr. and Robert W. Fuller. Mathematics of Classical and Quantum Physics. Dover, New York, 1992. [130] Supriyo Datta. Quantum Transport: Atom to Transistor. Cambridge University Press, New York, 2005. [131] Zhiben Ren. Nanoscale MOSFETs: Physics, simulation and design. PhD thesis, Purdue University, West Lafayette, U.S.A., October 2001. [132] R. Lake, G. Klimeck, R. C. Bowen, D. Jovanovic, D. Blanks, and M. Swaminathan. Quantum transport with band-structure and Schottky contacts. Phys. Stat. Sol. (b), 204:354-357, 1997. [133] J. W. Mintmire and C. T. White. Universal density of states for carbon nanotubes. Phys. Rev. Lett, 81(12):2506-2509, 1998. [134] Arkadi A. Odintsov. Schottky barriers in carbon nanotube heterojunctions. Phys. Rev. Lett., 85(1):150-153, 2000. [135] MPI - The message passing interface standard, http://www-unix.mcs.anl.gov/mpi, 2000. [136] Ramesh Venugopal, Zhibin Ren, and Mark S. Lundstrom. Simulating quantum transport in nanoscale MOSFETs: Ballistic hole transport, subband engineering and boundary conditions. IEEE Trans. Nanotechnol, 2(3):135-143, 2003. [137] Steven G. Louie. Electronic properties, junctions, and defects of carbon nanotubes. In Mil-dred S. Dresselhaus, Gene Dresselhaus, and Phaedon Avouris, editors, Carbon Nanotubes: Synthesis, Structure, Properties, and Applications, volume 80 of Topics in Applied Physics, pages 113-143. Springer, Berlin, 2001. [138] Enzo Ungersboeck, Mahdi Pourfath, Hans Kosina, Andreas Gehring, Byoung-Ho Cheong, Wan-Jun Park, and Siegfried Selberherr. Optimization of single-gate carbon-nanotube field-effect transistors. IEEE Trans. Nanotechnol, 4(5):533-538, 2005. [139] M. Pourfath, H. Kosina, B. H. Cheong, W. J. Park, and S. Selberherr. Improving DC and AC characteristics of ohmic contact carbon nanotube field effect transistors. In Gerard Ghibaudo, 104 Thomas Skotnicki, Sorin Cristoloveanu, and Michel Brillouet, editors, Proc. European Solid-State Device Research Conference, pages 541-544, Grenoble, France, September 2005. [140] T. Kamimura, K. Sakamoto, M. Maeda, K. Kurachi, and K. Matsumoto. Fabrication and I-V characterization of carbon nanotube single electron transistor operated at room temperature. In Proc. TMS 61st Annual Device Research Conference, pages 177-178, Salt Lake City, U.S.A., June 2003. [141] Jia-Lin Zhu, Hai-Feng Song, and Xiao Hu. Transverse-field and defect-azimuth effects in achiral carbon nanotubes. J. Phys.: Condens. Matter, 17:4629-4636, 2005. [142] L. Chico and R. Perez-Alvarez. Continuum model for long-wavelength optical phonons in cylinders: Application to carbon nanotubes. Phys. Rev. B, 69:035419-1-035419-6, 2004. [143] Gerald D. Mahan. Many-Particle Physics. Kluwer Academic/Plenum Publishers, New York, 3rd edition, 2000. [144] Mark Bockrath, David H. Cobden, Jia Lu, Andrew G. Rinzler, Richard E. Smalley, Leon Balents, and Paul L. McEuen. Luttinger-liquid behaviour in carbon nanotubes. Nature, 397:598-601, 1999. [145] A. Kasumov, M. Kociak, M. Ferrier, R. Deblock, S. Gueron, B. Reulet, I. Khodos, O. Stephan, and H. Bouchiat. Quantum transport through carbon nanotubes: Proximity-induced and intrinsic superconductivity. Phys. Rev. B, 68:214521-1-214521-16, 2003. [146] P. J. Burke. Liittinger liquid theory as a model of the gigahertz electrical properties of carbon nanotubes. IEEE Trans. Nanotechnol, 1(3):129-144, 2002. [147] S. Bellucci, J. Gonzalez, and P. Onorato. Doping- and size-dependent suppression of tunneling in carbon nanotubes. Phys. Rev. B, 69:085404-1-085404-6, 2004. [148] Calin Buia, Alper Buldum, and Jian Ping Lu. Quantum interference effects in electronic transport through nanotube contacts. Phys. Rev. B, 67:113409-1-113409-4, 2003. [149] Gordana Dukovic, Feng Wang, Dauhua Song, Matthew Y. Sfeir, Tony F. Heinz, and Louis E. Brus. Structural dependence of excitonic optical transitions and band-gap energies in carbon nanotubes. Nano Lett, 5(11):2314-2318, 2005. [150] Murray R. Spiegel. Mathematical Handbook of Formulas and Tables. Schaum's Outline Series in Mathematics. McGraw Hill, Toronto, 1968. 106 Appendix A Solution to Laplace's Equation For a closed cylindrical geometry, one may calculate the analytical solution to Laplace's equation with constant potentials applied to the various surfaces of a cylindrical box. The results presented here are standard, see Ref. [82] for example, but we present the derivation here for completeness. We solve V 2 V = 0 (A.l) in cylindrical coordinates subject to V(0, <f>, z) is finite, V(RG,<t>,z) = VG, V(p,<l>,0) = 0, V(p,<t>,LC) = VD, V(p,<j),z) = V(p,<t> + 2ir,z), where p = RG corresponds to the gate electrode, z = 0 is the source, z = L G is the drain, and we have chosen our potential reference here such that Vs = 0. Note that this solution computes the potential inside the cylindrical box before the nanotube is inserted, i.e., there is no charge inside the box. Using superposition, we can combine this solution with that for the case of homogeneous boundary conditions with charge inside, as derived in Appendix B. The Laplace problem will be solved as a linear superposition of two solutions. The first solution will consider the case where VQ = 0, and the second will be for VG = 0. Since each solution independently satisfies Eq. (A.l), the total solution will, too. It may be seen, as well, that the superposition will satisfy the desired boundary conditions listed above. 107 A . l Potential Applied to the Curved Surface The problem to solve, in cylindrical coordinates, is d_y idy yfFy &_y_0 ( A 2 ) dp2 p dp p2 d(j>2 dz2 subject to the boundary conditions mentioned previously. Using separation of variables,we assume a solution of the form V(p,4>,z) = R(p)G(4,)Z(z), (A.3) where R(0) is finite, R(RG) = VG, Z(0) = Z(LC) = 0, (A.4) 6(0) = 6(0 + 27r). (A.5) Substitution of Eq. (A.3) into Eq. (A.2) yields R" R' 6" -Z" -R+TR + 7Q = — (A'6) and since the left-hand side is independent of z, and the right-hand side is independent of 0 and p, each side must be equal to a constant, K2. This gives an equation for Z(z) as Z" + K2Z = 0, subject to the boundary conditions given by Eq. (A.4). The general solution is Z(z) = Ci sm(nz) + C2 COS(K.Z). 27(0) = 0 implies that C2 = 0, and Z{Lc) = 0 implies that n = mn/Lc, me [1,2,..., oo). The left hand side of Eq. (A.6) can be separated further to give where p, is another constant of separation like K was in the previous separation step. The equation for 0(0) is similar to the one for Z(z). The general solution is 0(0) = C 3 sin(£t0) + C 4 cos(/u0), (A.8) 108 where Eq. (A.5) implies that p, = I, I 6 [0,1,..., oo). The equation to be solved for R(p) is p2R"+pR'-(K2p2 + l2)R = 0. Making the substitution p = up gives p2R" + pR' - {p2 + l2)R = 0, which is Bessel's modified differential equation of order I [150]. The general solution is R{P) = C5I1(KP)+C6K1(KP), where I; and K; are the modified Bessel functions of the first and second kind, respectively, of order I. Since we require R(Q) to be finite, CQ = 0. Combining the results obtained thus far, and combining the constants, gives oo oo / \ / \ where Bim and Cim can be determined, using orthogonality, to be f J o L c ^ s i n ( ^ ) cos(^) dzd4> Bi„ / o 2 7 r / o i c V G s in (^ ) s i n (Z0 )d 2 d0 t Io" * ( ^ f f* ) ^n 2 ( T F ) sin2(l<t>) dzd<f>' Note that, in our case, since VQ is not a function of <f>, only the I = 0 term will contribute. Using this condition, and performing the integrations, noting that C; m drops out, gives ^ ) = £ B „ . o ( ^ ) where 2VG(1 - ( - 1 D B„ 109 A.2 Potential Applied to One End Again, we wish to solve Eq. (A.2), but subject to our second set of boundary conditions where the potential is zero everywhere on the boundary except on one end. As before, we separate variables, but for simplicity, we slightly change the separation constant, so that we get the first separation step to look like BT_ R_ _ -Z^_ _ _ 2 R + pR + p 2 B ~ Z ~ K ' In this case, the equation to solve for Z(z) is Z" - K2Z = 0. The general solution is Z(z) = Ci sinh(«;z) + C2 cosh(«;z), where C 2 = 0 so that Z(0) =0. . The second separation step looks like Eq. (A.7) except with a different sign on the n 2 p 2 term. o R R 9 9 @ 9 P l l + P R + K P = ^ = » -The solution for @(4>) is identical to Eq. (A.8). The equation to be solved for R(p) is now p 2 R " + p R ' + (K2p2 + l2)R = 0. Making the same substitution as before, p = up, gives p2R" + pR' + (p2 + l2)R = 0, which is Bessel's differential equation of order I [150]. The general solution is R(P) = CSJI(KP) + C6Y1(KP), where J; and Y; are the Bessel functions of the first and second kind, respectively, of order I. Since we require R(0) to be finite, C§ = 0. We also require that R(RG) = 0 for this case. This implies that K l m ~ R G ~ ' 110 where xim is the mth solution to Jj(a;/TO) = 0. Combining these results gives V(p, 0, s) = £ sinh (Blm sin(^) + Clm cos^)), where orthogonality can once again be used to determine our new Bim and Cim according t IoRG sin(Z0)pdpd<A Blm = C Io° ( ^ ) J 2 ( ^ ) sin2(^ )pdpd</>' Io'IORG ( ^ ) c o s d P d ^ > ' m " / 0 f l G sinh ( s g ^ ) J 2 (age) cos 2(Z0)pdpd0' If Vb has no 0 dependence, only the I = 0 term contributes, thus Bim drops out, and J o R G ^ J o ( ^ ) p d p C m — / o f l G ^ h ( ^ ) j 2 ( ^ ) p d p ' where the integrals may be performed analytically for the case where Vb is a constant. If we let XQmP then the numerator becomes R%VD f l O m / pJo(p) dp. Jo ' ' ' O m This is a standard integral [150] with solution B?GVD T , \ If we now let P p-Rc' then the denominator becomes R% sinh (^X°ft^C^ J pJl(x0mp) dp-This is also a standard integral which simplifies, since Jo{xom) = 0, to BQ . , (XomLc \ i T/ , N X 2 111 or ^ s i n h ( ^ c ) j 2 ( , 0 m ) . The final solution is OO / \ / / \ F ( p , < M = ] C C m J o ( ^ ) s i n h m=l ^ ' where x o m S i n h ^ a ^ j J i ( x 0 m ) A . 3 Total Solut ion The total solution to Laplace's equation, subject to our initial boundary conditions, VM.) - I [ « . ( ^ ) - ( ^ ) + ( ^ ) ( = g i ) ] where 2V G(1 - (-1)" 1) Dm. = m T r l o ( ^ ) ' x 0 m s inh(^22gi^j 112 Appendix B Green's Function for a Cylindrical Box A useful technique in solving Poisson's equation is to use a Green's function. Here, we derive this function for the case of a closed, cylindrical box with homogeneous boundary conditions. This solution may be added to that presented in Appendix A in order to solve Poisson's equation for the closed geometry cylindrical CNFET. The basic idea is to solve V 2 y = . g ( P - " ) . In cylindrical coordinates, this is 82V | ldV | 1 d2V | 9 2 V = g [ p ( p , 0 , z ) - n ( p , 6 z ) ] dp2 p dp p2 d(j>2 dz2 e Using the Green's function formalism, we can solve a different equation and use it to solve Eq. (B.l). We are seeking a function such that V(p, <f>, z) = - q- J \p(p', </>', z') - n(p', cf>', z')} GP(p, </>, z; p', cf>', z')p' dp' d</>' dz', (B.2) where the integral is over the whole space, and Gp is the Green's function. The appropriate equation for GP [82] is d2GP ldGp 1 d2GP d2GP 1., , . r / J / n n . - W + -p^F + 7 2 i W + ^ = P 6 { P ~ P ) 5 { 4 > ~ M Z ~ Z L ( B - 3 ) Gp, then, represents the impulse response of the system to a point charge. The Green's function for the potential inside a grounded cylindrical box, 0 < p < RG and 0 < z < Lc, may be found by solving Eq. (B.3) under the following boundary conditions Gp(0,(j),z) is finite, 113 GP{RG,4>,Z)=0, GP(p, 0,O) = O, GP(p,4>,Lc) =0, G P ( p , 0 , z ) = C7p(p>0 + 27T, z ) . This section will demonstrate the technique of eigenfunction expansions [120] for this problem. We look for the eigenfunctions and eigenvalues of the differential operator, and then express the Green's function as a linear superposition of these functions. If we let £ represent our operator, then Eq. (B.3) may be written as C G P = ^ 6(p-p')S(4>-4>')6{z-z'). (B.4) Let the eigenfunctions of C be Mm with their associated eigenvalues, A m . Since the Mm functions are a complete set, GP = J2NmMm(p,<t>,z), (B.5) m where Nm are constants determined by the boundary conditions. Eq. (B.4) then becomes J2 NmCMm = -5{p - p')5{4> - <P')5(z - z'). (B.6) Since Mm are eigenfunctions of C, CMm = ^mMm-Since they are also orthogonal, we can multiply both sides of Eq. (B.6) by the complex conjugate of one of the eigenfunctions, Mm, and integrate. If the functions are chosen to be orthonormal, we get Mm(p',<t>',z') = Nm\m, which can be rewritten and substituted into Eq. (B.5) to give GP(pA,z;p'^,z')=± Mmip>A\Z')Mm(pA,Z)^ ( B 7 ) m=l A m where we assume that none of the eigenvalues vanish. The eigenfunctions of the differential operator considered' here are found by solving d2M 1 DM 1 d2M d2M . . J n ,„ oX 114 subject to the same boundary conditions as for Gp above. Using separation of variables, we assume a solution of the form . ' M(p,<j>,z) = R(p)Q(<fi)Z(z). . Substitution of this form into Eq. (B.8) leads to the following set of equations in the usual way: Z" + K2Z = 0, (B.9) e" + /x20 = O, (B.10) • p2R" + pR' + [p2(-K2 - X) - p2] R = 0, (B.ll) where K and p are separation constants. The solution to Eq. (B.9) is Z(z) = Ci sin(Kz) + C2 COS(KZ), subject to 27(0) = 27(Lc) = 0. This implies that C2 = 0 and K = jn/Lc, j € [1,2,... ,00). Eq. (B.10) has a similar solution in form, but is solved under the constraint that 0(0) = Q((f> + 2ir). This is satisfied for 6(0) = C 3 e"* where the separation constant, p, must be an integer, I. For Eq. (B.ll), it is useful to make the substitution p = p\J—K2 — A in order to change the form of the equation to p2R" +pR' + (p2 - p2)i? = 0, which we recognize from Appendix A as Bessel's differential equation of order p. The general solution, noting that p, = I from the 6 solution, is fl(p) = C5J ;(p)+C6Y;(p). Since R(0) must be finite, CQ = 0. The other boundary condition is R(RG) = 0 which implies that RGV-K2 - A = Xim, where xj m has been defined in Appendix A. 115 With some algebra, we get the eigenvalues for M as , _ (xim\2 ( j n \ 2 Combining the above solutions, we get the eigenfunctions as oo OO f \ / • \ Mm(P,<f>,z)=Y: E c ^ j ' (^ ) s i n (F ) ' j = l l=—oo where Cji are constants that may be chosen using the orthogonality of eigenfunctions. In this case, we wanted to have orthonormal eigenfunctions, so these constants are chosen such that The z and <f> integrals are trivial, and the p integral can be rearranged to become one of those listed in Ref. [150]. The p integral is ^ { [ j ; ( ^ ) ] 2 + ( i - ^ ) j 2 ( ^ m ) ) , where the second term vanishes by the definition of xim. Using Ref. [150] again, this integral becomes R'Q T 2 / \ Substituting all of this work into Eq. (B.7) produces 2 ^ ^ ^ e " ( ^ ' ) j < J< ( ^ ) s i n ( £ ) s i n ( # ) GI' = -^RJX E E j=ii=-oom=i J,2+1(a;(m) in agreement with problem 3.23 in Ref. [82]. Note that only the I = 0 term will survive if there is no azimuthal variation. 116 Appendix C Conformal Mapping of the Cylindrical Laplace Equation We show how Laplace's equation changes under a conformal mapping using the methods of com-plex analysis [121]. This derivation will then be employed to transform the azimuthally-invariant cylindrical Laplace equation into an inverse space [71,84] as is useful when modeling open boundary CNFETs. First, we will start with some definitions. If we define a function g(U) = u(x, y) + iv(x, y) where U = x + iy, then we discover that, for g to be differentiate at some point, UQ = XQ + iyo, we require that the Cauchy-Riemann equations hold at that point. These are 5-5- ™ To see this, we can take the definition of the derivative and let £, a complex number, approach from any direction in the complex plane. If, for example, we approach along the horizontal direction, it is easy to show that 9 (Wo) = -Q^(xo,y0)+\—(xQ,y0). If we approach along the vertical direction, we get ,.,. . .du. . dv, . 9 ( w 0 ) = -i-z-(xo,y0) + —(x0,yo). oy dy These two expressions must be equal for the derivative to exist, therefore Eqs. (Cl) and (C.2) are necessary conditions for the function to be differentiable. It turns out that if the first partial 117 derivatives of u and v are continuous and satisfy the Cauchy-Riemann equations at all points in the domain, then g is analytic in the domain, i.e., it is differentiable everywhere inside the domain. See page 59 of Ref. [121] for a proof of this. C. l Transformation of the Cartesian Laplace Equation Suppose now that we have Laplace's equation d2V d2V _ dx2 dy2 ' where V is a real-valued function. The source of V may be represented as the real and imaginary parts of an analytic function. Since the real and imaginary parts of an analytic function have continuous partial derivatives to any order, d_du = d_du dy dx dx dy implies that S2v=_dh)_ dy2 ~ dx2' 1 ' where the Cauchy-Riemann equations have been used to get Eq. (C.4) from the elementary calculus result of Eq. (C.3). A similar result can be derived for u. Since the real and imaginary parts of an analytic function satisfy Laplace's equation, we can define a complex potential, and solve in the complex plane. This allows us to use conformal mapping to solve our problem. We now show how Laplace's equation transforms under a conformal mapping. The first important point is that the mapping must be one-to-one. This means that each point in the original W-plane must map to a unique point in the mapped space, the T-plane. If T = T(U) represents the mapping from one space to the other, this condition is T(Ui) = T{fA2) if and only if U\ = U2- It will be shown that, if T(U) is an analytic function itself, then Laplace's equation transforms in a very useful way. Namely, we recover Laplace's equation exactly in either domain. We will consider the specific case of the nu) = I transformation as this is of current interest, but the results and steps may be easily used with any other analytic mapping. First, we will confirm that this is an analytic function by checking that it 118 satisfies the Cauchy-Riemann equations, and that it is continuous everywhere. Note that we will consider a domain in the W-plane where the origin is not included, so this technicality will not be considered. Essentially, we consider the transformation of the space outside some disk into the interior of a different disk. 1 1 x — iy — = = —= ~ = u + iv U x + iy x2 + y2 implies that du x2 + y2 — 2x2 x2 — y2 dx (X2 + , , 2 )2 ( l 2 + y 2 ) 2 ' dv _ - (x2 + y2) + 2y2 _ x2 - y2 dy (x2+y2)2 {x2+y2)2' du 2xy dv dy =~(x2+y2)2 = 'fa' and hence the first partial derivatives are continuous and satisfy the Cauchy-Riemann equations, so T(U) is analytic. Following the steps outlined in problem 2 of Section 7.1 in Ref. [121], we let V(u,v) = V(x(u,v),y(u,v)) (C.5) where u and v are the coordinates in the T-plane. The gradient of V is given by V V - (— — )^ - — + i — \dx' dy J dx dy' where the subscript indicates that the gradient is performed in the W-plane. Similarly, yTV - (— — ^ - — + i — \du' dv J du dv' Using the chain rule on Eq. (C.5) yields W_dVdx dV_dy_ du dx du dy du dy_dydx + dV_dy dv dx dv dy dv This gives a relationship between VTV and V ^ ^ . Using the Cauchy-Riemann equations, V y + i f — — + — — dx dv dy dv \ dx dv dy dv 119 Rearranging and factoring yields where we can, again, use the Cauchy-Riemann equations to get where the star denotes a complex conjugate. This motivates us to define the gradient operator in the T-plane as ( d U Y ( 8 . d \ The Laplacian of a function, g, is defined by V 2 g — V • Vg, therefore we can write V\r as V2T = dU dT \dx2 dy2 / where we recall that a dot product of a vector with itself in the complex plane involves taking the complex conjugate of one of the arguments and multiplying. Finally, we have If then 82V d2V_(&V 82V\ du2 dv2 \ dx2 dy2 J d 2 V d 2 V A r . . w+W= f{x'y)' dU d T (C.6) (d2v d2v \ du2 dv2 d T dU : Nf(x(u,v),y(u,v)) may be solved in the transformed space. If Nf(x,y) = 0, then we can solve d2V d2V _ Q du2 dv2 In other words, Laplace's equation is invariant under the transformation. C.2 Transformation of the Cylindrical Laplace Equation This section transforms a three-dimensional Laplace equation, with azimuthal symmetry, into the space defined by 120 The equation we wish to transform is 82V ldV d2V _Q dx2 y dy dy2 ' where we recall the second derivatives may be transformed using Eq. (C.6). It now remains to transform the second term in this equation. First, we recall that u u2 + v2' v dVdx Wdy dx du dy du' dVdx dVdy dx dv dy dv' 2 / 9 . 9 1 uz + vz for a T — 1/U transformation. For convenience, we reproduce our starting four equations: dy du dV dv du dv dx dy' du dv dy dx We can multiply Eqs. (C.7) and (C.8) in order to help us eliminate the dV/dx terms: dxdy dv du dxdy du dv (C.7) (C.8) dV dx dx dV dy dx dx du dv dy du dv' dV dx dx dV dy dx dx du dv dy dv du Subtracting these gives dx dV dx dV _dV ( dy dx dy dx dv du du dv dy \ du dv dv du which can be rewritten using the Cauchy-Riemann equations as dxdy dv du dxdV du dv dy dy dx du + dy du This implies that dy dy dx du + dydy + dxdy du du du dv where Cauchy-Riemann has, again, been used. Putting this all together, and dividing out common factors, gives the transformed equation as 1 2 2 \ f92V d2V\ n dV <tt {c^ + l^)-2udu- + u2-v2\ dy dv = 0. 121 Appendix D The W K B Approximation Here, we provide details on using the WKB approximation for solving the one-dimensional, effective-mass Schrodinger equation following the derivation in Ref. [95]. We wish to solve h2 8 2w 2 ^ 8 ? = ™ - * * * ' ( D - 1 } for the case where we have an incident electron, of energy E and wavefunction ip, on a potential energy barrier U(z). Note that mes is the effective mass which is assumed constant throughout the solution space. It will be assumed that E < U(z) for z0 < z < z\, E > U(z) when z < ZQ or z > Z\, and E = U(z) when z = zn or z = z\. If we define z • z = —, U) 2m e f f w2' where w = z\ — ZQ, then we can rewrite Eq. (D.l) as ^ = (U(z)-E)ip. ' (D.2) The key to the WKB method is assuming a solution of the form V>(f) = JC( f )exp( i^ ) , where /C(z) and 9{z) are complex functions, ( < 1, and u is a constant. We further assume that /C(z) can be expressed in a series in powers of £, i.e., K[z) = /Co (z) + e/Ci (z ) + (z) + ... Substituting these into Eq. (D.2) and rearranging yields #c„ - £1 _ 2 t/c (e2f + 2\il-'K-ze-z + i ^ - ' / c ^ + (E-U)IC = O, 122 where a subscripted coordinate denotes a partial derivative with respect to that variable. Note that the exponential term cancels out of this expression. We now have an equation that can be solved for asymptotically. D.l Approximate Solution At this point, we note that by setting t = 1/2, and substituting our power series for /C, we get a sequence of equations in powers of £. Comparing terms that do not depend on £ gives ->C0(62)2 + {E-U)lCo = 0, which is known as the eikonal equation if we cancel the /Co out of this expression. Solving yields 6{z) = ± j' y/E-U{z)dz, where the lower limit of integration may be chosen to suit the problem. Comparing terms that involve yields -/Ci (02-)2 + 2i/C02<92- + i/Co02Z- + [E - U) /Ci = 0. This is referred to as the transport equation. Since 6 is a solution to the eikonal equation, this expression simplifies to /Cog _ Ogz /Co 2#2 with the solution where C is a constant. /Co = C ' \E-U\ D.2 Error In order to get a measure of the error, we can look at the next term in the expansion. Comparing terms involving powers of £ 2 t gives Kou - /C2 {Ozf + 2i/Ci2-02- + i/Ci02-2- + {E-U)JC2= 0, 123 which again simplifies, due to the eikonal equation, giving K-Ozz + 2i/Cif#z + i/Ci#2z = 0 . If we let K.\{z) — iK,o(z)W(z), and substitute for /Co in terms of 6, we get W — ^ s s s _ ^ (Qzz)2 which implies that W = C + \ ^ 2 - \ f ^ t d z . 4 ( 0 z ) 2 8J (0-z(z))3 We have a good expansion if |£/Ci| <S |/Co|, i.e., , if <C 1. If we write this condition in terms of E — U, we observe that the expansion breaks down near the turning points. We need to perform a more careful analysis about these points. D.3 Turn ing Points and Match ing In regions where the WKB expansion is valid, we can write down the general form of the wavefunction subject to the conditions specified above. This is ' ( ^ V f t e ' ^ ^ / ^ - ' + fte-'/iv^F-*) , z < z 0 *(*)=• ( ^ ) i ( c 8 e ^ V ^ « + C 4 e - / i > / m ^ d * ) , *><*<*. (srrjj <?5e J j ' V « , z > zi where Ci through C 5 are constants that need to be determined, zo = ZQ/W, and zi = zi/w. As we have already mentioned, these results are not valid close to zo and z\. If we examine the region close to zo, and expand E — U(z) about this point, we get E-U{z)~-U'{z0){z-z0), where the primed notation denotes a total derivative, and we have assumed that U' exists. In order to examine this region closely, we define z — zo = £8° C where Po is a constant and £ is an expanded coordinate representation. Substitution gives ?-afhi>« - u'{zo)e°w=0, 124 where we note that U'(zo) > 0 from our problem definition, and we choose BQ = 1/3 to simplify the equation. If we further define a0 = [U'(z0)}H, we get which is Airy's equation. The general solution is ip(a0) = C6Ai(a0) + C7Bi(a0), where C& and CV are constants, and Ai() and Bi(-) are the Airy functions of the first and second kinds, respectively. As a 0 —• oo, we can look at the asymptotic behaviour of these functions. Using the forms presented in Ref. [96], we get that ip{a0) r exp - - Q 0 2 + r exp -a 0 2 . Similarly, as an —» — oo, we have exp(i§(-a 0)*) exp(- i | ( -a 0 ) i ) , i i t , . —)=- (Cee-^ + C 7e^) + — L - ' (c 6e^ + C7e^) . In order to have a good solution, the asymptotic solution must join smoothly to the solution that is valid away from the turning point. If we recall that U(z) — E ~ U'(ZQ)(Z — ZQ), then 1 1 U-E ~ [U>(zo)]l {loo' Similarly, Near ZQ, for z > ZQ, we can write our WKB solution as 1 [U'(z0)}U^al Matching coefficients yields 2 A / 2 a C3exp ( -a§ 1 + C4exp I -g"o C 7 [t/'(f0)]^A 125 c 4 c 6 [U'(z0)}H& 2 ^ If we consider z < SQ near ZQ, we note that 1 1 E-U [U'(zo)}HH-a0) and J Zn E-U{z)iB__ IU'(ZQ) J za -A- 2 , ^ iz0 V "» V £ which provides the approximate form for the wavefunction in this region as 1 Ciexp ( -if;(-an)0 + C2exp A 2 ; (-an)2 [U>(z0)}-°^(-a0)k* I If we equate coefficients once again, we get Ci _ % - T = -1= (Cee-* + C 7e?) . This concludes the matching that is performed near the first turning point. If we examine the region close to Z\, we get E-XJ(z)~\U'(z{)\(z-z{). Again, we perform a stretching transformation via z — z\ = where we choose B\ = 1/3 as before. If we let «i = |c/'(«i)|*c, the equation to solve is with the general solution where Cg and CQ are constants. Now, as —ai —» oo, we have C 8 Ai(-ai)+C 9 Bi(-ai) , ^ ~ lJt^)\ G X P ( - ^ ( - Q l ) l ) + 0 F ( ^ ( ^ ( _ a i ) 126 Proceeding as before, we note that 1 1 U-E | t f ' ( f l ) | * £ j ( _ a i ) Note that 'z 0 V s J z a V S Jzl so we can define which is independent of z. Furthermore, so our wavefunction in this region is 1 | t f ' ( f i ) | ' £ & ' ( - a i ) i Performing the matching yields C3exp ( v- +C4exp (-*>.+ C 3 „ Cs -e | J 7 ' ( f i ) | * $ A 2 V ^ ' C4 _,, C i 9 e = | t f ' ( z i ) | * £ & A ' Finally, we must match as — a.\ —• —00 where the asymptotic form gives us exp I As before, 1 | t ^ i ) | 5 £ W The function to match to becomes 2 ! exp I l - a f 1 , I c ^ 1 1 q | * 7 ' ( z i ) | * V 3 127 which implies that 6 £ f j 2V7T V J \U'{zi)\H-> C8eT = - C 9 e ~ T . The matching is now complete, and we can solve for all but one of the unknown constants. The result, in terms of Ci , is i d ( 4 e 2 - - l ) 2 4e2" + l ' C3 = C4 = -C5 = C6 (4e2" + 1) (1 - i)' 4iC1V2e2u "(4e2" + l)( l - i ) ' 4C1e" 4 e 2 z / + 1 ' SiCiV^e 2 " [U'(z0)]H^ (4e2" + 1) (1 - i)' 2CiV^i : [U'izo)}*^ (4e2" + 1) (1 - i) 4CiV^re' / \U'(Zl)\i^ (4e2" + 1) (1 - i)' C 7 = C6 = C9 = , i |[/'(zi)| 5£A(4e 2" + l ) ( l - i ) Finally, we can compute the transmission probability T. If U(z) is constant outside of the turning points, and if we define U(z)=[ U S ' ~ Z < Z 0 ^ UD , z > zi then the transmission probability is given by T _ IE-UD | C 5 | 2 \ E-US C ! | 2 ' _ lE-Un ( 4e" V ~ lE-VD -2v where the last line is valid if e2" 3> 1. 128 Appendix E Solving for the Newton Update In Chapter 5, we formulated a quadratic matrix equation for the self-energy matrix representing a semi-infinite lead. Here, we show how one may compute the Newton update at each step of the iterative procedure, and how one may employ exact line searches in order to improve the convergence rate. E . l Computing the Raw Newton Update We consider the equation, from Chapter 5, r(g.) + r h r g s + [T^ST - (EI - ?L)} e = 0, (E.l) where we seek e. We define [125] P = QIT*ZU S = Q1[^g3T-(EI-rL)}Z1, f = Q2Z2, n = Q2~g\SZ2, where V, S, T, and It are, in this case, upper triangular, and Qi and Zj are unitary. These matrices may all be found via QZ factorizations, i.e., V, S, Q\, and Zi are found from a QZ factorization of rt and T^gsr — (EI — TL), and T, 71, Q2, and Z2 are found from a QZ factorization of the identity matrix and g\r^. Premultiplying by Q\ and postmultiplying by Q\ allows Eq. (E.l) to be written in the equivalent form [125] Q^Zxz\eZ2Z\rgsQ\ + Qx [T%T - (EI - TL)] ZXZ\EZ2Z\Q\ = -Qir(gs)Ql 129 Making the definitions y = Z\eZ2, and T = — Q\r{gs)Q\ yields vyv) + syr* = T, where we now seek y . If subscripts to a matrix represent the row and column entry, we may write this in the equivalent form N N m=l 1 = 1 where the star denotes complex conjugation. Since V, S, TZ, and T result from a QZ factorization, many of the terms in each sum vanish. Eliminating these terms yields Af N Fi,j = YL ('Pi,™ym,l'R-jtl + Si,mym,lTjti) • m=i l—j If we use the notation of Ref. [126], where a single subscript to a matrix denotes a column of the matrix, we can also write this as N N Tm = vYJ +s Yl y^,r j = m j = m Note that this differs slightly from the expression in Ref. [126], where the QZ factorization resulted in upper quasi-triangular, as opposed to upper triangular, matrices. Rearranging yields N {nmimP + T*tms)ym = Fm- Y (n*mJP + T m d s ) y h j = m + l which lends itself to a recursive computation [126]. E . 2 U s i n g E x a c t L i n e S e a r c h e s The idea behind exact line searches is that, instead of performing the usual Newton update of = _j_ £.) where the superscript denotes the iteration number, we update according to -(i+i) _ -(i) _j_ ^ w j j e r e o < < 2 is a damping factor. The interval for $ is a well-known restriction, required in order to obtain convergence with a damped Newton method, and may be derived by considering an appropriate norm of the residual at a particular iteration. The choice of i9 is determined by using additional information from the original quadratic matrix equation. 130 We have that r(gs + de) = r(gs) + •dr^erg, + 0 [SgsT - (EI - TL)} e + flVWe. We define p($) = \\r(gs + fie)]]2?, where \\A\\p = y/tiace(A^A) is the Frobenius norm, and choose $ such that a fourth-order polynomial, is minimized at each step with 0 < 1? < 2. It may be shown [125] that p($) either has an inflection point at •§ = 2 or a minima in (0,2], so such a d may always be found. 131
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Simulation studies of carbon nanotube field-effect...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Simulation studies of carbon nanotube field-effect transistors John, David Llewellyn 2006
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Simulation studies of carbon nanotube field-effect transistors |
Creator |
John, David Llewellyn |
Publisher | University of British Columbia |
Date Issued | 2006 |
Description | Simulation studies of carbon nanotube field-effect transistors (CNFETs) are presented using models of increasing rigour and versatility that have been systematically developed. Firstly, it is demonstrated how one may compute the standard tight-binding band structure. From this foundation, a self-consistent solution for computing the equilibrium energy band diagram of devices with Schottky-barrier source and drain contacts is developed. While this does provide insight into the likely behaviour of CNFETs, a non-equilibrium model is required in order to predict the current-voltage relation. To this end, the effective-mass approximation is utilized, where a parabolic fit to the band structure is used in order to develop a Schrödinger-Poisson solver. This model is employed to predict both DC behaviour and switching times for CNFETs, and was one of the first models that captured quantum effects, such as tunneling [i.e. tunnelling] and resonance, in these devices. In addition, this model has been used in order to validate compact models that incorporated tunneling via the WKB approximation. A modified WKB derivation is provided in order to account for the non-zero reflection of carriers above a potential energy step. In order to allow for greater flexibility in the CNFET geometries, and to lift the effective-mass approximation, a non-equilibrium Green’s function method is finally developed, which uses an atomistic tight-binding Hamiltonian to model doped-contact, as opposed to Schottky-barrier-contact, devices. This approach benefits by being able to account for both inter- and intra-band tunneling, and by utilizing a quadratic matrix equation in order to improve the computation time for the required self-energy matrices. Within this technique, an expression for the local inter-atomic current is derived in order to provide more detailed information than the usual compact expression for the terminal current. With this final model, an investigation is presented into the effects of geometrical variations, contact thicknesses, and azimuthal variation in the surface potential of the nanotube. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2010-01-18 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0065538 |
URI | http://hdl.handle.net/2429/18554 |
Degree |
Doctor of Philosophy - PhD |
Program |
Electrical and Computer Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_2006-199947.pdf [ 7.82MB ]
- Metadata
- JSON: 831-1.0065538.json
- JSON-LD: 831-1.0065538-ld.json
- RDF/XML (Pretty): 831-1.0065538-rdf.xml
- RDF/JSON: 831-1.0065538-rdf.json
- Turtle: 831-1.0065538-turtle.txt
- N-Triples: 831-1.0065538-rdf-ntriples.txt
- Original Record: 831-1.0065538-source.json
- Full Text
- 831-1.0065538-fulltext.txt
- Citation
- 831-1.0065538.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0065538/manifest