Self-consistent Vlasov-Poisson analysis of carrier transport in vacuum-based thermionic/thermoelectronic devices by Amir Hossein Khoshaman M.A.Sc., Simon Fraser University, 2011 B.Sc., Sharif University of Technology, 2009 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in THE FACULTY OF GRADUATE AND POSTDOCTORAL STUDIES (Electrical and Computer Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) August 2016 © Amir Hossein Khoshaman, 2016 ii Abstract Thermionic conversion involves the direct conversion of heat, including light-induced heat, from a source such as solar energy, to electricity. The progress of thermionic converters has been limited by issues such as the space charge effect and lack of materials with desirable mechanical, electrical and thermal properties. Nanotechnology could help address some of the main challenges that thermionic converters encounter. However, existing models, which were developed for macroscopic converters, are not adequate for many aspects of nanostructured devices. The work presented in this thesis primarily advances a new model to partially address this void and study emergent thermionic devices. We demonstrate a self-consistent and iterative approach to the Vlasov-Poisson system that overcomes the inherent limitations of the traditional methods. This approach serves as the foundation for more advanced and yet crucial cases of the operation of thermionic converters in the presence of back-emission, grids in the inter-electrode region and low-pressure plasmas. We develop the physics of the device in the presence of grids and demonstrate that momentum gaps could arise in the phase space of the electrons; taking into account these gaps, which had not been noticed in the past, is key to designing efficient thermionic converters and we predict improvements of 3 orders of magnitude in current density using a properly designed grid. We also develop the physics of the device in the presence of low-pressure plasmas, which are prime candidates for reducing space charge. We show that the output power density of a thermionic converter can improve by a factor of ~ 10 using a modest plasma pressure of 500 Pa. On a different front, we have also improved the traditional analytical model and developed an approach to extract the internal device parameters such as emission area and workfunction based on a limited set of experimental output characteristics. These parameters are highly dependent on the operating conditions and ex-situ measurements are not applicable. Therefore, our approach allows for a more systematic study of the device and material properties, which is key to further the development of thermionic converters, in particular based on novel materials and nanostructures. iii Preface The contributions from this thesis have led to the following journal paper publications: 1. A. H. Khoshaman and A. Nojeh, “Classical momentum gap for electron transport in vacuum and consequences for space charge in thermionic converters with a grid electrode,” J. Vac. Sci. Technol. B, vol. 34, no. 4, p. 40610, Jul. 2016. 2. A. H. Khoshaman and A. Nojeh, “A self-consistent approach to the analysis of thermionic devices,” Journal of Applied Physics, vol. 119, no. 4, p. 044902, Jan. 2016. 3. A. H. Khoshaman, A. Koch, M. Chang, H. Fan, A. Nojeh, “Nanostructured thermionics for conversion of light to electricity: simultaneous extraction of device parameters,” IEEE Transactions on Nanotechnology, vol. 14, no. 4, pp. 624–632, Jul. 2015. 4. A. H. Khoshaman, H. Fan, A. Koch, G. Sawatzky, A. Nojeh, “Thermionics, Thermoelectrics, and Nanotechnology: New Possibilities for Old Ideas,” IEEE Nanotech. Mag., vol. 8, pp. 4-15, 2014. The results of this thesis have also led to the following conference papers and presentations: 1. A. H. Khoshaman, A. Koch, A. Nojeh, “Light Induced Thermionic Energy Conversion using In-situ Intercalated Carbon Nanotube Forests,” The 16th International Conference on the Science and Application of Nanotubes (NT), Nagoya, Japan (2015) 2. A. H. Khoshaman and A. Nojeh, “A Comprehensive Approach to the Analysis of Nano-thermionic Converters through Particle Tracing,” The 28th International Vacuum Nanoelectronics Conference (IVNC), Guangzhou, China (2015) 3. A. H. Khoshaman and A. Nojeh, “Self-consistent Approach to Analysis of Nanostructured-Thermionincs,” The 16th International Conference on the Science and Application of Nanotubes (NT), Nagoya, Japan (2015) 4. A. H. Khoshaman, M. Chang, H. Fan, A. Koch, A. Nojeh, “Characterization of the internal parameters of nanostructured light induced thermionic emission devices for energy conversion,” The 14th IEEE International Conference on Nanotechnology (IEEE Nano), Toronto, Canada (2014) 5. A. H. Khoshaman, H. D. E. Fan, A. T. Koch, N. H. Leung, A. Nojeh, “Localized light induced thermionic emission from intercalated carbon nanotube forests,” The 27th International Vacuum Nanoelectronics Conference (IVNC), Engelberg, Switzerland (2014) iv 6. A. H. Khoshaman, H. D. Fan, and A. Nojeh, “Improving the Efficiency of Light Induced Thermionic Energy Conversion of Carbon Nanotubes by Potassium Intercalation,” The 97th Canadian Chemistry Conference and Exhibition, Vancouver, Canada (2014) 7. A. H. Khoshaman, M. Chang and A. Nojeh, “Extraction of Multiple Parameters of a Light-activated Thermionic Cathode with a Single Type of Experiment,” The 26th International Vacuum Nanoelectronics Conference (IVNC), Roanoke, USA (2013) 8. A. H. Khoshaman, M. Chang and A. Nojeh, “Laser Induced Structural Damage to Multi-walled Carbon Nanotubes in a Controlled-Pressure Environment,” The 57th International Conference on Electron, Ion, and Photon Beam Technology and Nanofabrication (EPIBN), Nashville, USA (2013) A list of other Journal and Conference papers to which I have contributed but are not directly related to the results of this thesis is provided in Appendix C. All of the Chapters were written by the author with the guidance of Professor Alireza Nojeh. All of the experiments and models were conceived by the author with guidance from Professor Alireza Nojeh; the experiments were performed with assistance from Mr. Harrison Fan, Mr. Andrew Koch, Mr. Mike Chang, Mr. Nate Leung and Dr. Mehran Vahdani Moghaddam. v Table of contents Abstract ........................................................................................................................................................ ii Preface ......................................................................................................................................................... iii Table of contents ......................................................................................................................................... v List of figures ............................................................................................................................................. vii List of symbols and constants .................................................................................................................. xii List of abbreviations ................................................................................................................................ xiv Acknowledgements ................................................................................................................................... xv Dedication ................................................................................................................................................. xvi Chapter 1 – Introduction .......................................................................................................... 1 Motivation ..................................................................................................................................... 1 Conversion of heat to electricity ................................................................................................... 2 The basics ...................................................................................................................................... 3 Factors affecting conversion efficiency ........................................................................................ 7 Reducing the emitter workfunction ............................................................................................... 9 1.5.1. Alkali metal intercalation .......................................................................... 9 1.5.2. Other approaches .................................................................................... 10 Mitigating the space charge effect .............................................................................................. 10 1.6.1. Applied electric and magnetic fields ....................................................... 11 1.6.2. Negative electron affinity materials ........................................................ 11 1.6.3. Reducing the inter-electrode spacing ...................................................... 12 Light-induced thermionic emission ............................................................................................ 13 1.7.1. CNT-based LITECs ................................................................................ 14 Theoretical background............................................................................................................... 17 1.8.1. Vlasov equation ...................................................................................... 17 1.8.2. Poisson equation ..................................................................................... 21 Survey of theoretical methods of modelling TECs ..................................................................... 21 Research objectives and overview ............................................................................................ 23 1.10.1. Inducement............................................................................................ 23 1.10.2. Overview ............................................................................................... 24 Chapter 2- An analytical approach to the analysis of thermionic energy converters ..... 29 Theory and model ....................................................................................................................... 29 Space charge limited regime ....................................................................................................... 32 Theory of the space charge limited regime ................................................................................. 32 Current-voltage characteristics at the critical and saturation points ........................................... 37 Output current-voltage characteristics in the space charge limited regime ................................ 38 Non-uniform temperature profile at the emitter .......................................................................... 40 Summary ..................................................................................................................................... 42 Chapter 3- A self-consistent approach to the analysis of thermionic devices .................. 43 Theory and model ....................................................................................................................... 43 Self-consistent solution in the absence of back-emission ........................................................... 44 vi 3.2.1. Particle tracing approach to bypass the analytical solution to the Vlasov equation ................................................................................................................ 49 Self-consistent solution in the presence of back-emission .......................................................... 52 Convergence of the solutions to the Vlasov-Poisson system ...................................................... 60 Summary ..................................................................................................................................... 62 Chapter 4- Analysis of thermionic energy converters comprising a grid ........................ 63 Theory and model ....................................................................................................................... 63 Vlasov equation in the presence of a grid ................................................................................... 63 Poisson equation ......................................................................................................................... 75 Results and discussion ................................................................................................................ 78 Summary ..................................................................................................................................... 90 Chapter 5- Low-pressure plasma-enhanced behaviour of thermionic converters .......... 91 Theory and model ....................................................................................................................... 91 Results and discussion ................................................................................................................ 94 Summary ................................................................................................................................... 103 Chapter 6- Contributions and future work ...................................................................... 104 Contributions ............................................................................................................................. 104 Future work ............................................................................................................................... 106 References ................................................................................................................................................ 108 Appendix A: Light-induced thermionic converter experimental results and comparisons with theory ....................................................................................................................................................... 113 A.1. Experimental set-up ................................................................................................................. 113 A.2. Extracting the experimental parameters ................................................................................... 115 Appendix B: Light induced thermionic energy conversion using in-situ intercalated carbon nanotube forests ...................................................................................................................................... 119 B.1. Results and discussion .............................................................................................................. 120 Appendix C: Publications that are not directly related to the results of this thesis .......................... 123 vii List of figures Figure 1-1. A schematic of a vacuum TEC: E is the emitter electrode, C is the collector electrode, and d is the inter-electrode distance. Electrons inside the gap constitute the space charge cloud. The device uses thermal energy from the heat source to deliver electric power to the load. .................................................. 3 Figure 1-2: The energy diagram of a vacuum TEC: 𝜙𝐸 is the emitter workfunction, 𝜙𝐶 is the collector workfunction, 𝑉𝑜𝑢𝑡 is the output voltage of the TEC, and 𝑒 is the elementary charge. The red area represents the thermal population of electrons. 𝐸𝑣𝑎𝑐, 𝑆𝐶𝐿 is the motive inside the inter-electrode space in the space charge-limited (SCL) regime and 𝐸𝑣𝑎𝑐 is the motive when the space charge effect is eliminated. ..................................................................................................................................................... 8 Figure 1-3. Comparison between (a) a regular TEC, where the entire area of the emitter is heated and (b) a light induced thermionic energy converter (LITEC) from MWCNTs, where only a small portion of the emitter is heated. ......................................................................................................................................... 15 Figure 1-4. The simulated values of the efficiency as a function of temperature for CNT-based LITECs with different workfunction values. Higher efficiencies can be achieved at lower temperatures as the workfunction of the material is reduced. This model includes space charge effects but neglects back-emission from the collector to the emitter. These curves were generated based on the model that was developed in Chapter 2 of this thesis. ......................................................................................................... 17 Figure 1-5. (a) Schematic diagram of a TEC. The device comprises for a hot emitter, and a collector, separated through a vacuum gap. The electrodes are connected externally through a load (represented by a bulb). (b) The motive diagram and the corresponding velocity distributions in the space charge limited mode of operation where a peak occurs in the motive profile (point 3). .................................................... 20 Figure 1-6. Schematic diagram of a TEC comprising a hot emitter, a collector, with vacuum in-between. The electrodes are connected externally through a load. A grid is inserted in the inter-electrode region and typically a positive voltage is applied to it with respect to the emitter. ...................................................... 24 Figure 2-1. (a) A typical I-V curve with the corresponding motive diagrams in each operation regime. (b)-(c) Typical motive diagrams corresponding to different operation regimes of a vacuum thermionic converter under negative (b) and positive biases (c). The subscripts “E” and “C” denote emitter and collector, respectively. 𝜓, 𝜇 and 𝜙 denote the motive, electrochemical potentials and workfunctions, respectively. ................................................................................................................................................ 30 Figure 2-2.(a) The dimensionless motive (𝛾) and dimension distance (𝜁) plotted as the solution of the coupled Poisson and Vlasov equations. (b) The values of the dimensionless motive, 𝛾 are fitted using linear, quadratic and cubic polynomial in the range 0 < 𝜁 < 15, as tabulated in Ref. [14]. These fitted polymonials are used to extrapolate the values of 𝛾 for 𝜁 > 15. As observed, having a limited range for dimensionless motive leads to substantial errors in case that a limited range of the values of the dimensionless variables are used. ............................................................................................................... 35 Figure 2-3. Comparison of the estimated value of the critical point voltage, 𝑉𝐶, when calculated using the assumptions made by Hatsopoulos, 𝑉𝐶, 𝐻𝑎𝑡, and the actual value of 𝑉𝐶 obtained from iteratively solving the problem at different (a) inter-electrode distances and (b) emitter temperatures. These results are calculated for typical values of device parameters of 𝜙𝐸 = 4.6 𝑒𝑉, 𝜙𝐶 = 4.0 𝑒𝑉. The radius of the heat spot, 𝑟, is assumed to be 100 𝜇𝑚. In part (a), 𝑇𝐸 = 2000 𝐾 and in part (b), 𝑑 = 700 𝜇𝑚. ..................... 36 Figure 2-4. The entire I-V characteristics and the important boundary points (critical and saturation voltages and currents) were calculated for a wide range of device parameters. The effects of altering (a) the inter-electrode gap size and (b) the hot spot temperature on the I-V characteristics are included in here. The values of the currents and voltages at the boundary points (critical and saturation points, viii represented by subscripts “C” and “sat”, respectively), are tabulated in here. For both (a) and (b), 𝜙𝐸 =4.6 𝑒𝑉 and 𝜙𝐶 = 4.0 𝑒𝑉, 𝑟 = 100 𝜇𝑚. For part (a), 𝑇𝐸 = 2,000 𝐾 and for part (b), 𝑑 = 500 𝜇𝑚. ........ 39 Figure 2-5. (a) Temperature profile (assuming linear temperature distribution) on the hot spot of the CNT-based LIETC. (b) Contributions of different temperature contour rings to the total current across a hot spot with a radius, 𝑟 = 200 𝜇𝑚. ................................................................................................................. 41 Figure 3-1. Flowchart representation of the proposed self-consistent approach to calculate the output characteristics of TECs. Initially, the electron density of the electrons is assumed to be zero, and therefore only the Laplace equation is solved. In later iterations, the Poisson and Vlasov equations, or alternatively, Poisson and the particle tracing equations are self-consistently solved to reach convergence. .................. 45 Figure 3-2. The evolution of the motive vs. distance as the loop in Figure 3-1 is repeated. 𝑛 represents the iteration number. The 𝑛 = 1 case corresponds to the Laplace equation solution. This initial attempt underestimates the space charge and therefore leads to high electron density. The exaggerated density results in higher space charge effect (𝑛 = 6), which in turn reduces the electron density in subsequent iterations. The parameters of the TEC under study are 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 2,000 𝐾 and 𝑑 =10 𝑚𝑚. ....................................................................................................................................................... 49 Figure 3-3. Comparison between the particle tracing approach to calculate the electron density and the analytical Vlasov solution. The average electron density (a) and motive (b) for the final round of the iteration, when the convergence is reached, are plotted here. In the particle tracing approach, 107 particles are tracked in each iteration. This plot corresponds to the converged solution after 20 rounds of iteration. The electron density and motive plots from particle tracing are calculated based on their average steady-state values. (c) Comparison between the improved-and-extended-Langmuir (IEL) and the numerical solutions (i.e. using Vlasov integrals or particle tracing) for the current-voltage characteristics of the same TEC as in part a) with 𝜙𝐸 =4.5 𝑒𝑉, 𝑇𝐸 = 2000 𝐾 and 𝑑 = 1 𝑚𝑚 (The overall current is calculated for an emission surface area of 3.14 × 10 − 8 𝑚2). The IEL results were obtained based on the strategy outlined in [11]. The particle tracing current is the average steady state and matches the solutions obtained by applying an asymptotic expansion to the electron density. ..................................................... 51 Figure 3-4. (a) The contributions of the emitter and collector current densities, 𝐽𝐸 and 𝐽𝐶, to the total current density, 𝐽𝑡𝑜𝑡 of a TEC with 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 2,000 𝐾, 𝑇𝐶 = 1,900 𝐾 and 𝑑 =1 𝑐𝑚. (b) The evolution of the total current density as a function of the iteration number, 𝑛. (c) The converged motive as a function of position for an applied voltage of 1 𝑉. (d) The contributions of the electrons originating from the emitter and the collector to the overall electron density as a function of position. ....................................................................................................................................................... 56 Figure 3-5. The characteristics of a thermionic converter in the presence back-emission. The workfunctions of the emitter and the collector are both 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, the emitter’s temperature is, 𝑇𝐸 = 1500 𝐾, the collector’s temperature is 𝑇𝐶 = 1400 𝐾 and the inter-electrode distance is 𝑑 =100 𝜇𝑚. Evolutions of the (a) current density from the emitter, 𝐽𝐸, (b) motive, 𝜓 (c) electron density, 𝑛𝐸, 𝑒𝑥 from the emitter and (d) electron density from the collector, 𝑛𝐶, 𝑒𝑥, are plotted. In part (a), 𝛿 is the mixing ratio between motive in subsequent iterations of solving the Vlasov-Poisson system ................... 59 Figure 3-6. The maximum output power density of a TEC as a function of the inter-electrode distance. The device parameters are 𝜙𝐸 = 2.5 𝑒𝑉, 𝜙𝐶 = 2.0 𝑒𝑉, 𝑇𝐸 = 1,800 𝐾 and 𝑇𝐶 = 700 𝐾. The proposed method was employed to calculate the output current density at different voltages. These values were subsequently used to calculate the power density. The maximum power density was found by sweeping the voltage and calculating the current at each value of the gap size, due to the space charge effect. The values calculated in this figure are in accordance with the values reported in the literature [14], [57]. ..... 60 Figure 4-1. (a) The potential (𝜓) or motive profile of a potential well as a function of distance, 𝑥 and the resulting velocity distribution functions, 𝑓𝑥, 𝑣𝑥 in different regions. (b) Schematic of a vacuum device comprising of an emitter, a collector and an ancillary electrode (grid) that allows the transport of electrons ix through the classical potential well. The potential well depicted in (a) can be generated due to the external applied fields as well as the interaction between the electrons in the inter-electrode region...................... 64 Figure 4-2. (a) The motive diagram of a regular TEC in the space charge limited mode. The possible velocity distribution functions are plotted as well. The electrons can maintain different velocity distributions depending on which side of the 𝜓𝑀 they are positioned at. The possible motive potentials and their corresponding velocity distribution functions are depicted in (b)-(h). (b) The point of maximum motive coincides with the point right outside the emitter. In this case, the electrons can only move unidirectionally. (c) The point of maximum motive occurs at the collector, and the space charge is not strongly felt. The possible velocity distribution functions are plotted underneath the motive profile. d) Space charge is more intense in the emitter-grid region. Electrons with enough energy to circumvent the maximum motive can move only unidirectionally. (e) The point of maximum motive resides on the collector, though the space charge is still more severely felt in the emitter-grid region. (f) The space charge is more intense on the grid-collector side and the maximum motive lies in this region. (g) The space charge is tangible in both emitter-grid and grid-collector regions, with the motive height at its maximum peak on the grid-collector side. (h) Same as g, save that the maximum peak occurs on the emitter-grid side and the peak on the grid-collector side serves as a local extremum point. ...................... 72 Figure 4-3. (a) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.001 for 𝑛 ≤ 500 and 0.01 for 500 < 𝑛 < 1000 and 0.1 for 𝑛 > 1000. The workfunctions are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, the emitter’s temperature, 𝑇𝐸, is 1800 𝐾. The collector’s temperature is assumed to be low enough so that the back-emission can be neglected. The inter-electrode gap distance, 𝑑, is 100 𝜇𝑚 and the grid is situated at 𝑥𝐺 = 10 𝜇. The bias on the grid, 𝑉𝐺, is 2 V. The applied voltage between the collector and emitter is 0 𝑉. b)The convergence of the motive in the IER for different values of the iteration number, 𝑛. (c) Evolution of the electron density as a function of position for different iterations. ................................................................................................................................ 80 Figure 4-4. The motive profile and the electron density of a TEC with 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, 𝑇𝐸 = 2,200 𝐾, 𝑑 = 100 𝜇𝑚 , 𝑥𝐺 = 10 𝜇𝑚 and 𝑉 = 0 𝑉. Also shown are the potential profile and the electron density for the same TEC, but in the absence of a grid. As detailed in the main text, the overall effect is a total reduction in the maximum motive barrier........................................................................................... 81 Figure 4-5. (a) The converged motive profile of a TEC with 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 3,000 𝐾, 𝑑 =1 𝑚𝑚, 𝑥𝐺 = 500 𝜇𝑚 for different values of the grid bias. (b) the resulting electron density as a function of position for different gate biases. (c) the maximum motive, 𝜓𝑀 as a function of the grid bias. The computed data can be estimated by a double exponential fit, 𝜓𝑀 = 𝑎 𝑒𝑥𝑝𝑏 𝑉𝐺 + 𝑐 𝑒𝑥𝑝𝑑 𝑉𝐺, where 𝑎 = 0.03798, 𝑏 = −3.983 , 𝑐 = 1.075 and 𝑑 = −0.2363 with 𝑅𝑀𝑆𝐸 = 0.01705. ...................................... 82 Figure 4-6. The maximum motive (red curves) and the current density (blue curves) as a function of the grid bias for different values of 𝑑. The current densities in the absence of the grid are plotted as straight lines. The parameters of the TEC are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, 𝑇𝐸 = 1,800 𝐾 and 𝑥𝐺 = 𝑑/2. ...................... 83 Figure 4-7 (a) The convergence of the motive in the inter-electrode region for different values of the iteration number, 𝑛. The workfunctions are 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, the emitter’s temperature, 𝑇𝐸, is 3,000 𝐾. The collector’s temperature is assumed to be low enough so that the back-emission can be neglected. The inter-electrode gap distance, 𝑑, is 1 𝑚𝑚 and the grid is situated half-way between the emitter and the collector, at 𝑥𝐺 = 500 𝜇. The bias on the grid, 𝑉𝐺 is 2 V. The applied bias between the collector and emitter is −2 𝑉. (b) Evolution of the electron density as a function of position for different iterations. c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.0001 for 𝑛 ≤ 8000 and 0.001 for 1200 > 𝑛 > 800 and 0.01 for 𝑛 > 1200. . 85 Figure 4-8. (a) The convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛, for a TEC with similar parameters as that of Figure 4-7. The applied bias between the collector and emitter is 0 𝑉. (b) Evolution of the electron density as a function of position for different x iterations. (c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.0001 for 𝑛 ≤ 8000 and 0.001 for 1200 > 𝑛 > 800 and 0.01 for 𝑛 > 1200. . 87 Figure 4-9. (a) The convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛, for a TEC with similar parameters as that of Figure 4-7. The applied bias between the collector and emitter is 2 𝑉. (b) Evolution of the electron density as a function of position for different iterations. (c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.0001 for 𝑛 ≤ 8000 and 0.001 for 1200 > 𝑛 > 800 and 0.01 for 𝑛 > 1200. ................................................................................................................................................... 88 Figure 4-10. (a) The convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛, for a TEC with similar parameters as that of Figure 4-9. The applied bias between the collector and emitter is −2 𝑉. The grid is positions in 𝑥𝐺 = 910𝑑 and 𝑉𝐺 = 10 𝑉 (b) Evolution of the electron density as a function of position for different iterations. (c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.001 for 𝑛 ≤ 1,000 and 0.01 for 𝑛 > 1000. ..................................................................................................................................... 89 Figure 5-1. The motive profile and the resulting distribution functions for electrons in two special cases where the number of valleys and peaks are less than 3 (the motive profile for positive ions is the negative of that of electrons). 𝜓 is the motive; 𝑓𝑥, 𝑣𝑥 is the velocity distribution function. 𝑣𝑥 represents the velocity along the direction of propagation (𝑥). Point 1 corresponds to the point just outside the emitter and the last point (12 or 11) represents the point just outside the collector. In the left motive, 3 regions with different velocity distribution functions exist, and there are two regions that have two different momentum gaps. For the motive on the right, there are 4 distinct regions and 3 different momentum gaps. These motives are due to the combined effects of both electrons and ions. For the velocity distributions at the top, the local minima of motive to obtain velocities 𝑣2 and 𝑣3 are given by the motive at points 1 and 4, respectively, whereas for the velocity distributions at the bottom, the local minima of motive to obtain velocities 𝑣2, 𝑣3 and 𝑣4 are given by the motive at points 1, 4 and 7, respectively. ................................. 93 Figure 5-2. The characteristics of a TEC in the presence of ions 𝛽 = 0, 10 − 4 and 5 × 10 − 4 . The workfunctions of the emitter and the collector are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, the emitter’s temperature is 𝑇𝐸 =1500 𝐾 and the inter-electrode distance is 𝑑 = 100 𝜇𝑚. The resulting (a) motive,𝜓𝑥 (b) electron density, 𝑛𝑒𝑥, and (c) ion density, 𝑛𝑖𝑥, are plotted. The case with 𝛽 = 0 represents the absence of ions. ............... 97 Figure 5-3. The characteristics of a TEC with similar parameters as that of Figure 5-2 but in the presence of ions with the flux density of ions being 1 × 10 − 3 of that of electrons. These parameters lead to oscillations in the output characteristics of this TEC. Evolutions (in iteration number) of the (a) current density from the emitter, 𝐽𝐸, (b) motive, 𝜓𝑥 (c) electron density, 𝑛𝑒𝑥, and (d) ion density, 𝑛𝑖𝑥, are plotted. The values of 𝛼 are depicted on the right axis in (a). .................................................................. 100 Figure 5-4. The motive and carrier densities of a TEC with 𝜙𝐸 = 2.5 𝑒𝑉, 𝜙𝐶 = 2.0 𝑒𝑉, 𝑇𝐸 =1500 𝐾, 𝑇𝐶 = 1300 𝐾, and 𝑑 = 100 𝜇𝑚. The ratio of the ion to electron flux density is set equal to 10 − 4 in both collector and emitter. The (a) motive and the total species and the total absolute density 𝑛𝑡𝑜𝑡 ≡ 𝑛𝐸, 𝑒𝑥 + 𝑛𝐶, 𝑒𝑥 − 𝑛𝐸, 𝑖𝑥 − 𝑛𝐶, 𝑖𝑥, (b) electron density resulting from the emitter, 𝑛𝐸, 𝑒𝑥 and collector, 𝑛𝐶, 𝑒𝑥 and (c) ion density resulting from the emitter, 𝑛𝐸, 𝑖𝑥 and collector, 𝑛𝐶, 𝑖𝑥, as a function of distance are plotted. ................................................................................................................ 102 Figure 5-5. The current density-voltage (𝑉𝑎𝑝) characteristics of the TEC with the parameters specified in Figure 5-2 under various voltage differences between the emitter and the collector. 𝛽 = 0 corresponds to the case where no ions are present. 𝜓𝑀 represents the maximum motive and 𝐽𝑡𝑜𝑡 is the total current density from the emitter to the collector. .................................................................................................. 103 Figure A-1. (a) The schematic of the experimental set-up for I-V measurements of a CNT-based LITEC. (b) The incandescent spot on the side-wall of the CNT forest as imaged by a CCD camera. (c) The xi internal parameters of the LIETC device are obtained by solving the reverse problem and fitting to the experimental results as depicted in the legend of the plot. ....................................................................... 115 Figure B-1. The sample holder apparatus used for the in-situ intercalation of carbon nanotubes under LITE conditions. Potassium is placed inside a custom-built electrically insulating crucible which is mounted on a resistive heater. The potassium is evaporated and deposited on the surface and side-walls of multiwalled carbon nanotubes. The entire set-up is placed inside a custom-built vacuum chamber. The output current-voltage characteristics are measured using a Keithley electrometer. ................................ 119 Figure B-2. The current-voltage characteristic of the MWCNT sample before the intercalation. The value of the workfunction of the CNTs was calculated based on our previous work [2] to be 4.6 𝑒𝑉. The input laser power for this experiment is 25 𝑚𝑊. .............................................................................................. 121 Figure B-3. The time evolution of current between the collector and the emitter in an intercalated sample under ex-situ conditions. ........................................................................................................................... 121 Figure B-4. The current voltage characteristics of the in-situ intercalated MWCNT forest under LITEC conditions. The incident laser power is 35 𝑚𝑊. The values of the workfunction of the emitter, the inter-electrode distance and the emitter temperature were calculated to be 3.5 𝑒𝑉, 980 𝜇𝑚 and 𝑇𝐸 = 1530 𝐾, respectively. .............................................................................................................................................. 122 xii List of symbols and constants General Table of Symbols and Constants Variable Symbol Description A Cross-sectional area 𝐴𝑅 Richardson constant 𝛿 Mixing ratio of motives 𝛽 Ratio of the saturation current density of ions to that of electrons 𝑑 Inter-electrode distance 𝑒 The unit charge 𝜖0 Permittivity of vacuum 𝑓(𝑟, 𝑣) Velocity distribution function 𝐼 current 𝐽𝐶 or 𝐽𝐶,𝑒 Current density of electrons from collector 𝐽𝑠𝑎𝑡,𝐶 or 𝐽𝑠𝑎𝑡,𝐶,𝑒 Saturation current density of electrons from collector 𝐽𝐸 or 𝐽𝐸,𝑒 Current density of electrons from emitter 𝐽𝑠𝑎𝑡,𝐸,𝑖 Saturation current density of ions from emitter 𝐽𝐸,𝑖 Current density of ions from emitter 𝐽𝑠𝑎𝑡,𝐸 or 𝐽𝑠𝑎𝑡,𝐸,𝑒 Saturation current density of electrons from emitter 𝐽𝑡𝑜𝑡 Total current density 𝑘𝐵 Boltzmann constant 𝐽𝐶,𝑖 Current density of ions from collector 𝐽𝑠𝑎𝑡,𝐶,𝑖 Saturation current density of ions from collector 𝑚𝑒 or 𝑚0 Mass of electrons 𝑚𝑖 Mass of ions 𝑛𝐶 or 𝑛𝐶,𝑒 The electron density from the collector 𝑛𝑚,𝐶 or 𝑛𝑚,𝐶,𝑒 The electron density from the collector at the point of maximum motive 𝑛𝐸 or 𝑛𝐸,𝑒 The electron density from the emitter 𝑛𝑚,𝐸 or 𝑛𝑚,𝐸,𝑒 The electron density from the emitter at the point of maximum motive 𝑛𝑚,𝐸,𝑖 The ion density from the emitter at the point of maximum motive 𝑛𝐸,𝑖 The ion density from the emitter 𝑛𝑚,𝐶,𝑖 The ion density from the collector at the point of maximum motive 𝑛𝐶,𝑖 The ion density from the collector 𝑛𝑡𝑜𝑡 The total charge density 𝜙𝐶 Workfunction of the collector 𝜙𝐸 Workfunction of the emitter 𝜓 motive 𝜓𝑀 or 𝜓𝑚 Maximum motive xiii 𝜓𝑀,𝑙 Local maximum motive 𝒓 Position vector 𝑟 Position 𝑇𝐶 Temperature of the collector 𝑇𝐸 Temperature of the emitter 𝒗 Velocity vector 𝑣 Speed 𝑣𝑥 x component of the velocity 𝑉 or 𝑉𝑎𝑝 Applied voltage between the collector and emitter 𝑉𝐺 The applied voltage between the grid and the emitter 𝑥 x coordinate of the position 𝑥𝑀 x coordinate of the point of maximum motive 𝑥𝑀,𝑙 x coordinate of the point of local maximum motive 𝑥𝐺 Position of the grid xiv List of abbreviations CCD Charged coupled device CNT Carbon nanotube GUI Graphical user interface HM Hemi-Maxwellian IEL Improved-and-extended-Langmuir ITO Indium tin oxide LITE Light induced thermionic emission LITEC Light induced thermionic converter MWCNT Multi-walled carbon nanotube NEA Negative electron affinity SCL Space charge limited TEC Thermionic energy converter UHV Ultra-high vacuum xv Acknowledgements I would like to express my heartfelt gratitude towards my supervisor, Prof. Alireza Nojeh, for his unwavering support, guidance, and advice throughout my PhD studies. He has been a wonderful mentor to me and has left an everlasting impact on me that will perennially morph how I will conduct my research in the future. One of his many laudable traits is that he furnishes his students with uttermost latitude and inexhaustible support, so that they could find the time to master the basics of their research first, rather than jumping headlong into a predefined path to produce immediate, albeit frivolous results. I would also like to thank Prof. Mona Berciu for introducing me to the astounding world of statistical mechanics and for her invaluable help and guidance during my studies. I am extremely grateful to Prof. George Sawatzky for fruitful discussions and his deep intuition. I am sincerely appreciative of Prof. Fabien Pease for his insights on my works. I am also thankful to my M.A.Sc supervisor, Dr. Behraad Bahreyni, and my B.Sc. supervisor, Dr. Bizhan Rashidian, for the myriads of things that I have learned under their supervision. I am very grateful to the Committee members for their time and also reviewing my thesis. I would like to thank my colleagues and friends, particularly Andrew Koch, Harrison Fan, Mike Chang, Mehran Vahdani Moghaddam, and Nate Leung for helping me with the experiments. I am genuinely thankful to my dear friend Rina Lai for helping me proofread parts of this thesis. xvi Dedication To my parents.1 Chapter 1 – Introduction1 Motivation Energy conversion determines many facets of modern civilization, including infrastructure, manufacturing, transportation, communication, etc. The average amount of available and harvestable energy per capita can be deemed as a telltale metric to gauge the progress of the human civilization. There is an ever growing demand for access to more efficient and higher density energy sources due to the rapid expansion of technology and the immutable proliferation of energy consumption. A testimony to this patent trend is the more than doubling of the global energy consumption in the last four decades [6]. This pattern is more tangible in the quickly expanding economies such as in China and India. For instance, China has seen a 300 % boost in energy consumption over the same time span [6]. This trend has induced irredeemable changes on the ecosystems and landscapes, which are exacerbated by emission of massive amount of carbon dioxide. This is primarily caused by burning fossil fuels and the exploitation of the planet’s resources. Moreover, there is a tremendous quantity of energy that is wasted in the form of heat. According to the United States Department of Energy, as much as 50% of the energy consumed by manufacturing processes is wasted as heat via different forms of heat transfer or hot exhaust gases and fluids [7]. Therefore, the development of more efficient and clean energy sources is indispensable to ameliorate the catastrophic effects that the fossil fuels have had on the planet. The significance of this problem and the potential for growth in this field can be easily appreciated by envisioning that the amount of the solar energy that floods the earth’s surface is so vast that, in about one year, it is twice larger than all the energy that can be harvested from all the Earth’s non-renewable resources (e.g., coal, oil, natural gas and even mined uranium) [8]. However, only about 1 % of the human energy sources are solar due to the expense of the development of solar harvesting technologies [6]. 1 Sections 1.2-1.7 of this chapter have been reproduced from [1] with minor modifications. Some parts of sections 1.8-1.9 are reproduced with minor modifications from [2] [3] [4] [5]. 2 Conversion of heat to electricity Harvesting the sun’s energy has the advantage of mitigating the negative effects associated with burning fossil fuels. At the same time, it is important to acknowledge the potential setbacks of using solar cells such as their intermittent power delivery and the ineluctable need for energy storage, the high cost of fabricating exotic materials that are usually key to high efficiency devices, the possible release of some of these potentially hazardous materials into the environment, and the need for vast swathes of land to be utilized to produce appreciable power. The term solar cell usually brings photovoltaics to mind. Indeed, the direct conversion of sunlight to electricity has dominated the realm of solar-energy harvesting for the past 50 years. However, researchers have now begun to re-evaluate the possibilities of thermionic and thermoelectric energy conversion because of some of their attractive features. These processes involve the conversion of heat to electricity, or light to heat to electricity, and have been investigated for over a century. Advances in nanoscale materials and fabrication techniques have now opened new doors for the application of these effects. The resurgence of interest in these areas is mainly because of two reasons. First, the ever-increasing amount of waste heat generated and the need for alternative, clean energy sources has provided a new impetus for seeking effective ways of capturing this abundant source of energy [9]. Second, advances in micro-/nanotechnologies have created opportunities for addressing some of the fundamental challenges that have plagued these devices since their inception [1] [10]. For example, the workfunction of the materials used in thermionic converters has a strong impact on the device efficiency. Nanomaterials present a new avenue for engineering structures with low workfunction that could still be stable at the high operating temperatures involved. Device dimensions and inter-electrode spacing are other critical factors. Here, again, micro-/nanotechnology enables previously unattainable dimensional control. Considering that sunlight is a rich source of clean energy, it is highly desirable to be able to harvest it using efficient heat engines that do not involve sophisticated mechanics and moving parts. Thermionic and thermoelectric generation are excellent candidates for this purpose. However, reaching the required operational temperature is far from trivial due to the issues related to heat spread in the conductive electrodes of traditional devices. Recently, it has been shown that carbon nanotubes (CNTs) can be used as electrodes that, while being electrically conductive, can effectively trap heat and attain very high temperatures with low incident optical intensities, thus creating a new path for light-activated thermionics [11]. Hence, some of the unique properties of nanostructures may offer new solutions to issues once thought insurmountable, such as the intricate relationship between electrical and thermal conductivity in metals. 3 Here, we review the basics of thermionic and thermoelectric conversion. The important factors involved in optimizing the performance of converters based on these mechanisms are then discussed. In each case, we describe examples of solutions offered by micro-/nanostructures. The basics A photovoltaic cell relies on the excitation of electrons through the absorption of photons and the subsequent spatial separation of electrons and holes. Fundamentally, the excitation and separation of electrons and holes does not necessarily have to be induced by light. For example, another way to excite electrons to higher energy levels is by heating the material. A thermionic energy converter (TEC) converts heat directly into electric power by means of thermionic electron emission, thus acting as an electric heat engine and without using moving parts like steam turbines. Typically, TECs are comprised of two main electrodes, as depicted in Figure 1-1. It is noted that sometimes the acronym TIC, rather than TEC, is used in the literature to denote thermionic converters. Figure 1-1. A schematic of a vacuum TEC: E is the emitter electrode, C is the collector electrode, and d is the inter-electrode distance. Electrons inside the gap constitute the space charge cloud. The device uses thermal energy from the heat source to deliver electric power to the load. The electrons are thermionically emitted from the hot electrode (emitter or cathode) into a vacuum (or some other medium), traverse the inter-electrode distance because of their kinetic energy, and, eventually, are collected at the cold electrode (collector or anode). A negative charge thus builds up on the collector, gradually hindering further electron collection, until, eventually, the net electron flux between the electrodes becomes zero. If the external circuit is completed by connecting an electric load between the two electrodes, the charge gathered at the collector will flow back to the emitter through the load, generating a steady-state output voltage and current. 4 Thermoelectric conversion is, in essence, similar to thermionic conversion, but without a vacuum gap separating the emitter and collector. Although thermoelectric effects were first discovered in metals, modern-day thermoelectric devices consist of semiconducting p-n junctions. In the Seebeck effect, both electrical contacts to the load are cooled, and the junction is kept at a higher temperature than the contact ends (or vice versa). If the temperature gradient is increased, the charge carriers can more readily surmount the junction potential barrier. The electrons travel from the n- to p-layer, and the holes from the p- to n-layer, resulting in an electromotive force. Alternatively, when a current passes from one material to the other, the kinetic energy of the electrons is changed, and the difference appears as heating or cooling at the junction; this is the related, but opposite, Peltier effect. The terms potential and motive are used interchangeably in this thesis. This term was introduced by Langmuir and is a scalar whose negative gradient at any point is equal to the force exerted on an electron at that point. In the cases studied in this thesis, this corresponds to the negative of the electrostatic potential generated by a steady-state charge density. Thermionic conversion may be viewed as a thermodynamic steam engine cycle that uses electrons as the working fluid. In a TEC, the emitter can be thought of as the electron boiler, while the collector is the electron condenser, leading to an electrical pressure (potential) gradient across a load. This will produce work similar to that done by vapor pressure in steam engines. In the ideal case, the overall system efficiency approaches the Carnot efficiency, since electron evaporation may be considered a reversible process, and the temperatures of the hot and cold electrodes may be assumed to be constant during the process. Additionally, irreversible mechanisms such as friction are smaller compared to the situation in engines working with fluids, due to the absence of macroscopically moving parts. The main sources of irreversibility are the transport of electrons in an electric field (a non-equilibrium process), the Joule heating inside the load and the thermalization of the electrons at the anode. The thermionic current can be described by the Richardson–Dushman equation. This equation is found by summing up the contributions of all electrons having a velocity component normal to the emission surface and outward, assuming Fermi–Dirac statistics. By introducing several approximations, which are valid in most practical applications, the overall result is 𝐽 = 𝐴∗𝑇2 exp (−𝜙𝑘𝐵𝑇) , (1-1) where 𝐽 is the current density, 𝑇 is the absolute temperature, 𝜙 is the workfunction of the material, 𝑘𝐵 is the Boltzmann constant, and 𝐴∗ is the apparent emission constant of the material [12]. A number of early landmark efforts laid down the foundations of TEC technology. Although the quantitative description of thermionic emission was first presented by Richardson in 1902, the discovery of 5 this phenomenon by Edison can be traced back to as early as 1885. Nonetheless, the earliest practical demonstrations of thermionic emission for energy conversion purposes emerged in the mid-1950s, when breakthroughs occurred in the technologies of high-temperature-stability electrodes and powerful heat sources [13]. The inter-electrode gaps of about 100 μm in these TECs were formed by means of precision machining. These devices achieved efficiencies of 10–15% [14] and were subsequently used as power sources for space missions in the following two decades. However, the inherent electronic efficiency of a TEC, or the efficiency that strictly corresponds to the electronic processes under ideal electronic transport conditions, is much higher. In principle, this conversion efficiency can approach 90% of a reversible heat engine (the Carnot efficiency) [14], which can reach values of more than 60% under ideal conditions. Because of the exponential dependence of TEC current on temperature, small increases in temperature can lead to substantial improvements in efficiency. Thermionic converters have the potential to operate at exceedingly high temperatures—on par with temperatures generated by the burning of fossil fuels. On the other hand, turbines require operation at much lower temperatures to preserve structural integrity, due to the presence of mechanical stress and vigorous conditions such as hot fluids and chemical products resulting from the combustion process. Therefore, the maximum attainable efficiency is considerably higher in thermionic converters and closer to the Carnot efficiency. For instance, coal is one of the most dominant fuel sources—responsible for the generation of a substantial portion of the world's electricity—and burns at about 1,500 °C, whereas turbines usually operate significantly below this point (about 700 °C) [15]. TECs exhibit inherent advantages over their turboelectric counterparts due to the lack of moving parts in their structure. This feature endows them with potential advantages in cost, reliability, life-time and weight. Since the new surge of interest in nanomaterial-based TECs has only recently started burgeoning, the exact amounts of their cost and life-time are not well-estimated. However, due to the remarkable properties of nanomaterials and their inherently higher thermal and mechanical stabilities, it can be anticipated that these emergent devices’ performance would surpass that of their antecedents from 1960s and 1970s. Thermionic converters can have exceptionally high power densities compared to other technologies. Typical values of power densities on the order of 104 𝑊/𝑚2 have been reported since 1960s for TECs [16]. Higher values of power density are easily conceivable, considering the exponential dependence of the current density on temperature; for instance, for an emitter workfunction of 2.0 𝑒𝑉 and an emitter temperature of 2,000 𝐾, a power density of ~ 107 𝑊/𝑚2 can be achieved if the space-charge is completely mitigated. This simple example illustrates the momentous impact of reducing the space-charge effects and reducing the workfunctions of the electrodes to harness the remarkable potential of TECs. For comparison, the power density generated by the most efficient mass-produced photovoltaic solar cells is ~ 175 𝑊/𝑚2 [17]. This is partly due to the fact that a higher power density would require a larger incident photon 6 intensity which could have a catastrophic effect on the performance of the photovoltaic device due to increase in temperature. In contrast, the power density of a light-activated thermionic converter can only increase by increasing the intensity of the light. The power densities achieved by windmills are slightly higher than in photovoltaics and can have values as high as 400 𝑊/𝑚2 [17]; for instance, for an industrial 70 MW turbine, the swept area is 23 × 104 𝑚2 [18] . Plantation of fast-growing hybrid poplars deliver power densities that are 100 times lower than this and are typically below 1 𝑊/𝑚2 [9]. Gas turbines, when combined with steam turbines, exhibit the highest power densities among these competing technologies. This is due to high efficiency of electricity generation and higher densities of fuel extraction, making power densities as high as 103 𝑊/𝑚2 realizable [9]. A comparison between TECs, thermoelectric generators and steam turbines is provided in Table 1-1. TEC Thermoelectric generator Steam turbine Efficiency of current devices ~ 15 − 20 [19] ~ 6 − 7 % [20] ~ 50 − 55 % [21] [22] Projected efficiency ~ 0.90 𝜂𝐶 a Unclearb ~ 50 − 55 % Initial cost High High High Maintenance cost Low Low High Size Compact Compact Bulky Energy carriers Electrons and ions Electrons, holes and phonons Water vapour Moving parts No No Yes Limiting factor Space-charge; workfunction; thermal stability; heat dissipation ZT; thermal stability Maintenance cost; large infrastructure; not able to operate at high temperatures; highly irreversible Temperature dependence Exponentially increasing current density with temperature Material degradation and a smooth increase in efficiency Cannot operate at high temperatures a Carnot’s efficiency b The efficiency of thermoelectric generators could ideally approach Carnot’s efficiency for infinitely large ZT. However, ZT>3 has never been achieved [23]. Table 1-1. Comparison between some technologies that involve direct conversion of heat to electricity. 7 Factors affecting conversion efficiency The conversion efficiency of TECs can be defined as the ratio of the power density produced in the load to the input flux of heat on the emitter, 𝑞𝐻, where the subscript H stands for hot, in accordance to the notations employed in physical chemistry texts [24]. The power per unit area in the load is equal to (𝑉 − 𝑉𝐿𝑒𝑎𝑑)𝐽, where 𝑉 is the difference between the electrochemical potentials of the emitter and collector, 𝑉𝐿𝑒𝑎𝑑 is the voltage across the leads, and 𝐽 is the total current density. There are excellent sources in the literature that discuss the efficiency of thermionic converters and limits on how far their efficiencies can be improved [25] [26] [27]. As will be explained in section 1.10, this thesis is primarily concerned with filling gaps in the existing models for the transport of charge carriers in TECs to account for effects that could become important in modern devices based on nanomaterials and nanotechnology. Calculation of conversion efficiency would involve the knowledge of the specifics of a complete converter device structure and assumptions for loss through the leads, heat conduction to surroundings and other issues that are not of immediate relevance to the transport of charge carriers in vacuum. Describing those effects, which have already been studied in detail in the literature, would constitute a major detour from the scope of this thesis. In this section, however, for the sake of completeness, we provide a brief overview of the loss mechanisms in TECs (from their cathode) and refer the reader to the references cited above for more detail. 𝑞𝐻 can be considered as the sum of several terms. The first one is heat flux due to the transport and emission/collection of electrons: 𝑒𝐽 (𝜙𝐸 + 𝜓𝑀) + 2𝑒 𝑘𝐵(𝑇𝐸 𝐽𝐸→𝐶 − 𝑇𝐶 𝐽𝐶→𝐸), where 𝐽𝐸→𝐶 and 𝐽𝐶→𝐸 determine the current from the emitter to collector and vice versa, and 𝜓𝑀 is the maximum motive. It is noted that we have not included a term 𝑒 𝐽𝐶→𝐸 𝜙𝐶, in accordance with the thermodynamics definition of the efficiency in a Carnot cycle [24]. The next term is due to heat loss flux through conduction to the surroundings and the leads. The last term is due to heat loss flux through thermal radiation of the emitter (in the far-field regime equal to 𝜎𝜖𝐸(𝑇𝐸4 − 𝑇04), where 𝜎 is the Stefan-Boltzmann constant and 𝜖𝐸 is the emissivity of the emitter). In a more intricate calculation, the thermally radiated power reabsorbed by the emitter from the collector should also be included. The heat loss through radiation is an abruptly increasing function of temperature, making low-temperature performance of TECs more favorable. Nonetheless, due to the exponential dependence of current density on temperature, the efficiency becomes very large at extremely high temperatures (though the material stability would be a major issue). The two most significant factors affecting the efficiency of TECs are the workfunctions of the electrodes and the space charge effect in the inter-electrode region. Conceptually, the workfunction of the emitter can be considered the energy barrier to the thermal evaporation of electrons into free space, as seen in the TEC energy diagram in Figure 1-2. According to the Richardson–Dushman equation, the number of electrons that are able to surmount this barrier is proportional to exp(−𝜙𝐸/𝑘𝐵𝑇𝐸) for a given electrode size 8 and duration of time, where 𝑇𝐸 is the temperature of the emitter. The difference between the emitter and collector workfunctions, 𝜙𝐸 − 𝜙𝐶 can be considered a measure of the open-circuit voltage. Figure 1-2: The energy diagram of a vacuum TEC: 𝜙𝐸 is the emitter workfunction, 𝜙𝐶 is the collector workfunction, 𝑉𝑜𝑢𝑡 is the output voltage of the TEC, and 𝑒 is the elementary charge. The red area represents the thermal population of electrons. 𝐸𝑣𝑎𝑐,𝑆𝐶𝐿 is the motive inside the inter-electrode space in the space charge-limited (SCL) regime and 𝐸𝑣𝑎𝑐 is the motive when the space charge effect is eliminated. Ideally, the collector workfunction should be as small as possible, and, as a rule of thumb, the emitter workfunction should be about 1 eV higher than the collector workfunction to maximize efficiency [14]. So far, no suitable metals have been found that combine stability at high temperatures and sufficiently low workfunctions for near-optimal performance. To maintain a stable, low workfunction, the electrodes are often immersed in alkali metal vapor, due to the low workfunction of alkali metals. The electrons ejected from the hot emitter have a finite speed and thus require a finite amount of time to travel to the collector. The electrons occupy the inter-electrode space during this time, forming a cloud of negative space charge. This gives rise to an electric field, which repels newly ejected electrons, and only those with sufficient kinetic energy to overcome this repulsion may reach the collector. This effect can be seen in the energy diagrams in Figure 1-2. In vapor TECs, this electric field is reduced, to a degree, by the positive ions of the alkali metal present in the vapor. However, some electrons collide with these vapor 9 atoms elastically and inelastically, which may cause them to return to the emitter, adding to the undesirable reverse current [14]. Reducing the emitter workfunction The emitter materials of TECs have traditionally been chosen from refractory metals with workfunctions on the order of 4–5 eV. This leads to the need for very elevated temperatures to achieve reasonable levels of emission current and output power. These high temperatures lead to implementation difficulties and increased losses due to the spread of heat to the surroundings and incandescence from the hot emitter. Therefore, the metals were typically coated with low-workfunction materials, which were not stable. 1.5.1. Alkali metal intercalation One way to circumvent this issue in nanomaterials is through the inclusion of alkali metals into the small spaces between atomic/molecular layers of the host, also known as intercalation. This process requires energy to increase the distance between the host layers against the van der Waals force. This energy can be supplied by means of external heat as well as exothermic interaction due to charge transfer between the guest molecules and the host. Two common examples include stage-1 (C8K) or stage-2 (C24K) K/CNT, where the stage number refers to the number of potassium atoms between the carbon layers, C is carbon, K is potassium, and CNT stands for carbon nanotube. Intercalated nanomaterial emitters have potential advantages over refractory metal emitters with alkali metals deposited on the surface. The latter are relatively unstable and tend to continuously evaporate or change their morphology on the surface, as they have low melting and boiling points. Therefore, a reservoir containing alkali metals is usually necessary to have the TEC operate consistently [13]. Intercalated nanomaterials, on the other hand, not only substantially reduce the workfunction of the host material but may also exhibit higher temperature stabilities [28]. Various types of carbon have been studied as thermionic emitters [29]. Therefore, developing methods of reducing their workfunctions is highly desirable. Pristine carbon nanofibers and nanotubes have workfunctions in the range of 4–5 eV, similar to that of polycrystalline graphite [30]. The workfunction of a CNT depends on several factors, such as its chirality, diameter, and surface oxidation condition [31]. However, by optimization of these parameters, the workfunction can be only marginally reduced. It is also possible to modify the workfunction of CNTs using adsorbates. However, because of weak intermolecular 10 attraction, the stability of the resulting structure is not adequate for the high temperatures required by TECs. On the other hand, it is possible to significantly reduce the workfunction of CNTs by alkali metal intercalation. It has been observed that intercalates of potassium or cesium with single-walled CNTs exhibit workfunctions of 3.3 and 2.4 eV, respectively [32] [15]. Intercalates have also been formed from stoichiometric reactions of molten potassium with various types of graphitic carbon nanofibers. Stage-1 K with herringbone graphitic carbon nanofiber (GCNF) showed a reduced workfunction of about 2.2 eV and stability up to 1,300 K [31]. Other work has shown that potassium intercalated CNTs can be stable up to 820 K [15]. 1.5.2. Other approaches Another interesting method to reduce the workfunction is by means of nanostructure geometry engineering. By coating the surface of the emitter with a thin metal layer with periodic ridges, a series of quantum wells are formed [33]. By imposing the resulting additional boundary conditions on the electron wave function at the surface, it is possible to introduce forbidden states. Electrons are banned from filling these states and, therefore, occupy higher energy levels, raising the chemical potential and correspondingly reducing the workfunction. Advances in stable low-workfunction materials have prompted investigations into their use as emitters in vacuum TECs. N-type diamond thin films, doped with nitrogen, have exhibited effective workfunctions of less than 2 eV. Low-electron-affinity aluminum gallium nitride (AlGaN) thin films also show promise as TEC emitter materials, although band bending due to surface states has pushed the workfunction to above 2.3 eV in practice [34]. The energy barrier that electrons see before thermionic emission can be effectively reduced in a gas discharge environment. The local electric field at the emitter can be modified because of the presence of ions exiting the discharge region. Go et al. theoretically studied the influence of ions from a gas discharge on thermionic emission of electrons. Their analysis reveals that the presence of ions can significantly increase the thermionic emission current, with more tangible effects in the case of stationary ion charges [35]. Mitigating the space charge effect As mentioned previously, the cloud of electrons between the two electrodes creates a potential (motive) barrier for the electrons ejected from the emitter. This effect is exacerbated as the current density increases, as the distance between the electrodes increases, and as the average electron speed decreases [14]. Several solutions have been proposed to mitigate the space charge effect. 11 1.6.1. Applied electric and magnetic fields One way to decrease or eliminate the potential barrier created by the negative space charge cloud is to apply an electric or magnetic field to the electrode gap via an auxiliary electrode. This effectively allows one to engineer the electron distribution in space to minimize the space charge barrier. The first of the two most-prominent triode TEC designs is the magnetic triode, in which crossed electric and magnetic fields steer electrons toward a collector lying on the same plane as the emitter. The other configuration is the electrostatic triode, in which a grid is placed between the electrodes to change the electric potential landscape to accelerate electrons. In vapor TECs, the negative space charge cloud is meant to be compensated for with positive ions. These ions can be produced by thermionic emission from the surface of the emitter or by collision of emitted electrons with vapor atoms. However, the ionization potential of these alkali metals is usually around 4 eV, and most of the thermionically emitted electrons do not have sufficient kinetic energy to ionize these atoms. Therefore, the ions are usually introduced by some other means, such as an arc discharge between the emitter and an auxiliary electrode, where the electric field is strongest. This results in a significant loss in the output power of the TEC. More recently, advanced vapor TECs have been proposed that use an auxiliary grid electrode in combination with a longitudinal magnetic field. These TECs theoretically reduce or eliminate the space charge effect. Two device configurations have been proposed by Moyzhes et al. [36]. In the first one, the hot electrons needed for ionizing the alkali metals are trapped in a potential well and separated from the thermionic current. In the second one, the alkali metal vapor atoms are ionized directly on the gate electrode. 1.6.2. Negative electron affinity materials It is possible to lower the electrons' potential energy level just outside the emitter surface—the vacuum level of the emitter—to fall below its conduction band. This class of emitters is known as negative-electron-affinity (NEA) materials. Several surface orientations of diamond exhibit NEA behaviour when terminated by hydrogen. This property, along with the low thermionic barrier in doped diamond, makes it an attractive candidate for thermionic emission. This mechanism opens a route to alleviate the space charge effects. When electrons emitted from the conduction band reach the vacuum level, their kinetic energy will noticeably increase, since the vacuum level rests at a lower energy than the conduction band. Therefore, NEA can effectively reduce the number of slower electrons, which spend the most time in the inter-electrode region and thus contribute most significantly to the space charge barrier [30]. NEA materials may also present an effectively lower barrier to electron emission due to quantum-mechanical tunneling. 12 1.6.3. Reducing the inter-electrode spacing Decreasing the distance between the electrodes is a straightforward approach to mitigate the space charge effect, although energy loss via radiative heat transfer between the electrodes may become significant, especially as the electrode spacing approaches the submicrometer range. At large gaps, heat is transferred between the electrodes primarily via propagating electromagnetic waves, following the Stefan–Boltzmann law. Evanescent electromagnetic waves, which decay exponentially away from the electrode surface, i.e., surface plasmon polaritons, are also excited by incident radiation. When the electrode gap shrinks to a distance on the order of the characteristic wavelength of blackbody radiation, these waves begin to couple to each other in an effect sometimes called photon tunneling [37]. In TECs with this small of an electrode gap, this form of radiative heat transfer becomes rapidly dominant and severely limits the efficiency. Lee et al. [38] have combined the equations for SCL thermionic current with far-field and near-field heat-transfer calculations to determine the optimal emitter–collector gap for a TEC. For near-ideal emitter and collector materials, the optimal gap for peak efficiency was found to lie between 0.5 and 1 μm, depending on emitter temperature. At electrode gaps of this size or smaller, the space charge effect essentially becomes insignificant. Fortunately, devices with such electrode gaps can be readily fabricated using current microelectromechanical systems fabrication technologies. Thermionic devices with gaps as small as 1.7 μm have been fabricated with silicon carbide (SiC) emitters and silicon collectors using standard surface micromachining techniques [10]. Thermionic conversion was successfully demonstrated by resistively heating the suspended SiC electrode up to 2,000 K (in the 50-μm-gap device), and no undesired contact was made between the electrodes. However, the conversion efficiency was still severely limited by the high workfunction of SiC. The previous discussion concerning the efficiency of TECs with submicron electrode gaps only dealt with thermionic emission of electrons. In theory, another peak in efficiency could potentially occur at an electrode gap on the order of a few nanometers, when electron tunneling becomes significant. The characteristics of nanoscale-gap TECs based on cesium iodide-coated graphite have been studied based on electron tunneling and radiative-heat-transfer calculations [39]. The electric current was found to be drastically higher but with a comparably severe increase in evanescent-wave radiative heat transfer. For an electrode gap of this size, the electrode materials would need to be carefully chosen for their optical properties to minimize evanescent-wave coupling [40]. Another significant challenge lies in the fabrication of devices with such small gaps. Electrodes with a vacuum gap smaller than 5 nm and area larger than 7 mm2 have been successfully fabricated using an 13 electroplated copper/silver/titanium/silicon (Cu/Ag/Ti/Si) structure, where the layers separate as a result of differing thermal expansion coefficients and controlled adhesion properties [41]. A different proposed solution is to place dielectric nanowires between planar electrodes to obtain a fixed, nanometer-scale gap [42]. Some of the highest efficiencies of TECs have been achieved by using high temperature high pressure plasma enhanced TECs with refractory metals such as molybdenum, tantalum or tungsten at high cesium pressures. Calculated efficiencies on the order of ~ 20 % have been reported for molybdenum electrodes at temperatures up to ~ 1980 𝐾 [19]. Light-induced thermionic emission The energy of photons can be exploited to heat the emitter material to sufficiently high temperatures and thus have light-induced thermionic emission (LITE). Naturally, an abundant and free source of light is the sun. A solar thermionic generator is a special case of a TEC, namely, a light-induced TEC (LITEC). In a single-junction photovoltaic cell, photon absorption is limited to the bandgap energy, which creates an inherent limitation on efficiency. Energy from incident photons exceeding the bandgap will be lost as heat, while photons with sub-bandgap energies are simply not absorbed, severely restricting the quantum efficiency. This effect is mitigated in multi-junction photovoltaic cells, which absorb photons at several discrete energies. A TEC relies on heating, which does not directly depend on this bandgap energy, and, theoretically, the entire absorbed portion of the spectrum of the incident light maybe used. Moreover, the structure of a LITEC is, in principle, relatively simple and robust and does not require high-quality semiconductors. LITECs are not well established in the clean-energy industry yet due to several practical issues in their implementation. The intricate connection between thermal and electrical conductivity is primarily responsible. An excellent electrical conductivity is necessary in a TEC, but this normally comes with high thermal conductivity, which leads to substantial heat loss to the surroundings. Consequently, extremely high incident powers and large, elaborate light-focusing contraptions are required to reach the optical intensities required to achieve the desired temperature. Heating the emitter to temperatures of >1,500 K remains challenging, even with the most advanced technology. Of note is NASA's High Power Advanced Low Mass solar thermionic system, which features extremely large, complex light-collection and -focusing structures that are often several meters in diameter [43]. 14 In light of the difficulties associated with solar heating of the emitter, most tests performed on LITEC devices have used alternative methods such as electrical heating. Under carefully designed experimental conditions, efficiencies as high as 11% have been achieved in this manner [43]. Light-harvesting TECs convert solar radiation into thermal energy, producing electricity via a heat engine, whereas photovoltaic cells use the quantum nature of light to excite electrons to higher energy bands. The processes could potentially be combined through thermally enhanced photoemission or photon-enhanced thermionic emission [32] [44] [45] [46]. The situation may further be improved by coating the surface with metal nanoparticles, which has been shown to enhance light absorption via a plasmonic process by increasing the effective optical path length inside the active layer, therefore increasing overall absorption [47]. 1.7.1. CNT-based LITECs Nanomaterials exhibit intriguing thermal and mechanical properties and could potentially circumvent some of the challenges associated with traditional TECs. Several salient features, such as high surface area, high mechanical strength, and resilience to high temperatures, make them promising candidates for thermionic applications. CNTs are particularly interesting in this context because of additional desirable characteristics, such as their quasi-one-dimensionality and high absorptivity over a broad spectral range. Electrons are confined in the transverse directions, limiting the energy states into which they can scatter, thus increasing the mean free path and conductivity. 15 Figure 1-3. Comparison between (a) a regular TEC, where the entire area of the emitter is heated and (b) a light induced thermionic energy converter (LITEC) from MWCNTs, where only a small portion of the emitter is heated. 16 Interestingly, CNTs have recently been shown to overcome the fundamental challenge of the spread of light-induced heat in emitters. When an array of multiwalled CNTs (a so-called CNT forest) is illuminated by a sufficiently focused low-power beam of light, a heat-trap effect is observed. This highly localized heating mechanism allows the illuminated spot to be heated to thermionic emission temperatures (>2,000 K), without significantly heating the surroundings [11]. This effect was attributed to two main factors, including a rapid drop in the thermal conductivity versus temperature in CNTs, as opposed to a less-rapid drop in bulk materials, as well as the quasi-one-dimensional nature of heat transport in CNTs as contrasted with isotropic bulk materials. By mitigating the heat spread to the surroundings due to the localized nature of the heat spot, this effect eliminates a major loss mechanism in thermionic conversion. Among other parameters such as the density and uniformity of the CNT forest, the threshold incident power required to induce the heat-trap condition depends primarily on the area of the illuminated spot. For an incident light beam diameter of several hundred micrometers, the threshold intensity is on the order of 50 W/cm2 [11]. Assuming that photons of different energies of the solar spectrum all contribute to heating, the average sunlight intensity on the surface of the earth, once focused with a handheld glass lens, will suffice to thermionically emit electrons from CNTs based on this effect. Accordingly, a solar LITEC has been demonstrated based on this phenomenon [48] [49]. The efficiency of the first device at peak power was low (10−6%). However, this low efficiency is not a fundamental limit (as shown by NASA's LITEC devices [43]), and major improvements are possible. Moreover, even this early prototype exhibited a short-circuit current density and peak power density comparable to those of state-of-the-art photovoltaics. This is a testament to some of the inherent advantages of LITECs, which are understood to enable higher power densities than photovoltaic devices. Moreover, the open-circuit voltage of these devices can also be appreciably higher than that of most photovoltaic devices. The most conspicuous improvements will arise from a reduction of the workfunction of the CNTs. Although the issue of heat transfer to the surroundings is minimized, major energy loss takes place through incandescence from the hot spot. By reducing the workfunction of the CNTs and thus operating at lower temperatures (without sacrificing the electron emission current), the energy loss due to radiation of heat can be substantially reduced. As the theoretical prediction of Figure 1-4 shows, the efficiency can be enhanced by several orders of magnitude if the workfunction of CNTs is reduced by about 2 eV. The results of Figure 1-4 were generated based on the model that was developed under the works of this thesis (Chapter 2). 17 Figure 1-4. The simulated values of the efficiency as a function of temperature for CNT-based LITECs with different workfunction values. Higher efficiencies can be achieved at lower temperatures as the workfunction of the material is reduced. This model includes space charge effects but neglects back-emission from the collector to the emitter. These curves were generated based on the model that was developed in Chapter 2 of this thesis. Further improvements can be attained by reducing the space charge effect through the optimal design of the device as well as by using aspherical lenses and more-sophisticated optics to reduce chromatic aberrations and allow for the better focusing of sunlight. Theoretical background The major part of this thesis comprise of solving Poisson and Vlasov equations self-consistently. In this section we provide an overview of these equations and how they can be employed in the context of thermionic converters. 1.8.1. Vlasov equation An understanding of the Vlasov equation can directly follow from a discussion on the Boltzmann transport equation. The classical theory of transport processes is described by the Boltzmann transport equation. The solution to this equation results in a distribution function in the phase space, 𝑓(𝒓, 𝒗). The 18 distribution function is a function of six variables of the phase space (𝒓 is the position vector and 𝒗 is the momentum vector). The number of particles in a differential volume (𝑑𝒓𝑑𝒗) in the phase space (𝒓, 𝒗) is equal to 𝑓(𝒓, 𝒗) 𝑑𝒓𝑑𝒗. The Boltzmann transport equation follows directly from the Liouville theorem of classical mechanics, which asserts that the distribution function remains constant when a volume element is followed along a flow-line [50] [51] [52]. Therefore, we have 𝑓(𝑡, 𝒓, 𝒗) = 𝑓(𝑡 + Δ𝑡, 𝒓 + Δ𝒓, 𝒗 + Δ𝒗), where t represents time and Δt represents a change in time. In the presence of collisions, this equation is modified to: Δ𝑡 (𝜕𝑓(𝑡, 𝒓, 𝒗)𝜕𝑡) |𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 = 𝑓(𝑡 + Δ𝑡, 𝒓 + Δ𝒓, 𝒗 + Δ𝒗) − 𝑓(𝑡, 𝒓, 𝒗), where (𝜕𝑓(𝑡,𝒓,𝒗)𝜕𝑡) |𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 represents the rate of change of the number of particles in the element with respect to time due to collisions [53]. Now we can consider the definition of derivative (as Δt → 0) and use the chain rule to obtain the Boltzmann transport equation: (𝜕𝑓(𝑡, 𝒓, 𝒗)𝜕𝑡) |𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 =𝜕𝑓(𝑡, 𝒓, 𝒗)𝜕𝑡+ 𝒗. 𝛁𝒓𝑓(𝑡, 𝒓, 𝒗) +𝑭(𝑡, 𝒓, 𝒗)𝑚. 𝛁𝒗𝑓(𝑡, 𝒓, 𝒗). In the steady state, we have (𝜕𝑓( 𝒓, 𝒗)𝜕𝑡) |𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 = 𝒗. 𝛁𝒓𝑓( 𝒓, 𝒗) +𝑭( 𝒓, 𝒗)𝑚. 𝛁𝒗𝑓( 𝒓, 𝒗). The time-independent Vlasov equation can be obtained from setting the collisional term equal to zero in the Boltzmann equation: 𝒗. 𝛁𝒓𝑓( 𝒓, 𝒗) +𝑭′( 𝒓, 𝒗)𝑚. 𝛁𝒗𝑓( 𝒓, 𝒗) = 0. (1-2) The central insight is that the external forces on the particles 𝑭( 𝒓, 𝒗) are changed to 𝑭′( 𝒓, 𝒗) to include the inter-particle forces as well as the external forces. This equation is not self-sufficient and is coupled with the Poisson equation (section 1.8.2) to evaluate the forces as the gradient of the motive in a mean-field fashion. In a 1-D case and in the absence of magnetic fields, the steady-state collisionless Vlasov equation can be written as 19 𝑣𝑥𝜕𝑓(𝑥, 𝑣𝑥)𝜕𝑥−1𝑚𝑑 𝜓(𝑥)𝑑 𝑥𝜕𝑓(𝑥, 𝑣𝑥)𝜕𝑣𝑥= 0, (1-3) where 𝑥 denotes the direction normal to the emitter’s surface, 𝜓(𝑥) represents electron motive (in this case, electric potential energy), 𝑣𝑥 is the velocity of the electron along the direction of propagation and 𝑓(𝑥, 𝑣𝑥) is the velocity distribution function in one dimension [18] [54]. Note that the 𝑦 and 𝑧 coordinates are cyclical [55] and the dependence of the distribution function on their corresponding velocities is trivial. The general solution to this equation can be written as 𝑓(𝑥, 𝑣𝑥) = ∑𝑓𝑖 (𝜓 +12𝑚 𝑣2) ,∞𝑖=1 where 𝑓𝑖 is any arbitrary function of 𝜓 +12𝑚 𝑣2. 𝑣 is the speed of the electrons. Langmuir [56] proposed to write this solution in the form 𝑓(𝑥, 𝑣𝑥) = ∑𝑓𝑖 (𝜓 +12𝑚 𝑣2)∞𝑖=1= ∑𝑎𝑖 exp {−𝑏𝑖 (𝜓 +12𝑚 𝑣2)}∞𝑖=1. (1-4) The constants 𝑎𝑖 and 𝑏𝑖 can be calculated by imposing the correct boundary conditions at a certain point in the phase-space of the electrons. So far, the discussion of the 1-D Vlasov equation has been general. When electrons originate from one region and travel to another region (this is the case in TECs as portrayed in Figure 1-5), the correct boundary condition is that electrons adopt a hemi-Maxwellian distribution at the point of maximum motive (potential): 𝑓(𝑥𝑀, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32exp (−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥), (1-5) where Θ represents the unit step function and 𝑛(𝑥𝑀) is the electron density at 𝑥 = 𝑥𝑀 (point 3 in Figure 1-5. (b)). The motivation behind this assumption is that electrons originate from one electrode (emitter) and therefore those with kinetic energies lower than 𝜓𝑀 cannot surmount the maximum motive [57]. In section 3.2 we will provide arguments on the validity of this assumption in the case of the emergent TECs and also in the presence of back-emission. 20 Figure 1-5. (a) Schematic diagram of a TEC. The device comprises for a hot emitter, and a collector, separated through a vacuum gap. The electrodes are connected externally through a load (represented by a bulb). (b) The motive diagram and the corresponding velocity distributions in the space charge limited mode of operation where a peak occurs in the motive profile (point 3). 21 In the following chapters of this thesis, this equation is solved under different schemes that could arise in a TEC (see section 1.9 for a brief overview). 1.8.2. Poisson equation The second important mathematical underpinning of this work is solving the Poisson equation ∇2𝜓(𝑟) = −𝑒2𝑛(𝑟)𝜖0, where 𝑛(𝑟) is the electron density at point 𝑟. Poisson equation is pertinent to this problem in the absence of time-varying magnetic fields (which is the case in typical TECs). Also, a question could arise to the validity of the static Poisson equation to the patently dynamic problem that we are studying. A Full treatment of a problem requires solving the entire Maxwell’s system of equations. However, within the context of thermionics, since the attained velocities of the particles are much smaller than the speed of light, this quasi-static approximation is completely valid [58] [59] [60]. In a 1-D case, Poisson equation can be written as 𝑑2𝜓(𝑥)𝑑𝑥2= −𝑒2𝑛(𝑥)𝜖0. (1-6) It is highlighted here that this system of equations can provide a good description of realistic experimental conditions. For instance, the theoretical calculations of this thesis agree with the analytical solutions obtained by Langmuir (in cases where Langmuir’s model can provide a solution), which have been shown to match well with the experimental results [14]. As another example, as shown in the Appendices, the results of this thesis can be used to study the behaviour of thermionic energy converters with MWCNTs as the emitter. Survey of theoretical methods of modelling TECs A method pervasively employed for the analysis of TECs was put forward several decades ago by Langmuir, Hatsopoulos and Gyftopoulos [14] [56]. It involves solving the Vlasov-Poisson system of equations in the space charge limited (SCL) regime, assuming the electrons in the inter-electrode space follow the statistics of a steady-state collisionless gas (in an external potential that involves their Coulomb 22 interactions) with a hemi-Maxwellian (HM) velocity distribution at the point of highest motive [14] [57]; this gives rise to a closed-form solution in the form of an integral, which can then be calculated numerically. This approach has a narrow range of applicability in the presence of strong space charge, and suffers from high splicing errors due to employment of different equations at different operation regimes. In this method, the solutions only exist when the motive reaches a maximum in the inter-electrode distance. It is because analytical solutions, in a closed form, arise when the boundary conditions are enforced such that the derivative of the motive is zero at one point in the inter-electrode region. This condition significantly narrows down the applicability of the analytical solution in the modes of operation that are not space charge limited. This effect is most palpable at the boundaries between the retarding mode and the space charge mode, where the highest errors can occur if a combination of precise numerical integration and an iterative strategy are not employed [2]. This method was promulgated by Hatsopoulos’s landmark book on thermionics and has been adapted by many researchers ever since [14] [57]. Smith et al. [61][62] developed the theory of TECs with a negative electron affinity (NEA) emitter using an approach similar to that of Langmuir and Hatsopoulos, and proposed a method to calculate the limits of the space charge regime. Smith also proposed that the space charge effects can be mitigated in the case of a TEC incorporating an NEA material as the anode [63]. Lee et al. have modeled the behaviour of TECs with very small gaps where near-field heat transfer becomes important and have observed that the optimum gap roughly corresponds to the characteristic wavelength of the emitter's thermal radiation [38]. Another approach was adopted by Stephanos and Meir et al. [15] [14], which involves solving the differential equation that is obtained by replacing the analytical solution to the Vlasov equation, inside the Poisson equation, and employing a non-linear solver. This method can treat all regimes of operation and overcome the splicing errors; however, the solution is limited to the cases where an analytical solution to the Vlasov equation exists, e.g., when the TEC comprises of two infinitely large and uniformly heated flat electrodes. A finite difference time-domain (FDTD) particle-in-cell (PIC) approach to obtaining the output characteristics of thermionic converters has been adopted by Lo et al. using a 1-dimensional object oriented particle tracing software package called OOPD1[65]. This method can be computationally expensive, if one is interested in the steady state behaviour of the system, as is usually the case in energy conversion applications, since the electrons are followed in real time and their influence on each other is accounted for by solving the Poisson equation. 23 Research objectives and overview 1.10.1. Inducement As touched upon in the previous sections, and as will be discussed further in here and the following sections, although experimental efforts are favorably underscored to build highly efficient thermionic converters, meager attention has been focused on the development of explanative models. The quintessence of this trend can be appreciated by considering this revealing example: Hatsopoulos proposed that the usage of a positively-biased grid electrode could reduce space charge in his classic text (Figure 1-6). Since then, this scheme has been discussed on many occasions in the literature, notably in recent times in works such as [36] and [15] and especially using sub-nanometer-thin, electron-transparent grids based on graphene [66]. Using this limited model misses the subtle but crucial phenomenon of the momentum gap that we have discovered. The approach and results we present in Chapter 4 solve this problem correctly and are indispensable in tackling the fundamental challenge of space charge and improving the efficiency of TECs. Indeed, our results demonstrate that by proper design of a grid, the current density in the device can be augmented by 3 orders of magnitude. 24 Figure 1-6. Schematic diagram of a TEC comprising a hot emitter, a collector, with vacuum in-between. The electrodes are connected externally through a load. A grid is inserted in the inter-electrode region and typically a positive voltage is applied to it with respect to the emitter. 1.10.2. Overview The goal of this thesis is to develop explanative models that can describe the operation of TECs under various operating conditions and regimes. The first step is improvement of the current analytical approach to expand its sphere of applicability to the wide range of parameters that one encounters in the emergent thermionic devices. The method proposed by Hatsopoulos and Langmuir can be used to evaluate the current-voltage (I-V) characteristics in the space charge limited mode of operation, assuming the device parameters are known. However, this strategy has never been used to solve the reverse problem of obtaining the internal parameters of the device from the I-V experimental curves. This is because conventional thermionic devices are made using metal electrodes with known properties, such as workfunction, and the entire electrode is typically heated; therefore, parameters such as emission spot area are well known and 25 temperature can be easily measured. Hatsopoulos and Gyftopoulos assumed a uniform temperature profile in their model [14][57]. However, as will be shown in Chapter 2 and appendix A, a large temperature gradient may exist in nanomaterials-based thermionics. For instance, in a CNT-based light induced thermionic energy converter (LITEC), only a small portion (with a radius, 𝑟, on the order of several hundred micrometers) of the nanotube array may be heated. The exact value of the workfunction depends on the local properties and temperature of the illuminated spot [67]. Moreover, the presence of adsorbates can have noticeable effects on the workfunction of CNTs. For instance, the adsorption of water molecules will reduce the workfunction [68] [69] while the presence of organic molecules on the CNT tip will result in larger workfunctions [67]. Additionally, the local density of states at the tip and sidewalls of CNTs are different, leading to different values of workfunctions [70]. Therefore, the study of these emerging thermionic devices requires a strategy to characterize their internal parameters. In addition, the lack of sufficient computational power and a stable algorithm have made it difficult to tackle the more computationally intensive reverse problem in the past. The integrals involved in this method do not have exact analytical solutions and behave non-linearly as the internal parameters of the device are varied. In fact, the numerical approximations tabulated in Ref. [14] are sufficient to predict the space charge limited behaviour of thermionic converters only in a limited range of temperatures, workfunctions, and inter-electrode spacings. Although we have improved the analytical method of analysis of TECs to solve the reverse problem and have adopted a self-consistent approach to minimize the splicing errors (Chapter 2), there is an inherent limit to the “reach” of the analytical method (in Deutsch’s terminology [71] [72]). In more elaborate cases, e.g. when back-emission of carriers from the collector is not negligible, or when a grid is inserted in the inter-electrode region, or when the TEC is operating in a low-pressure plasma, it is not possible to find analytical solutions to the Vlasov-Poisson system. Therefore, we propose a different strategy in Chapter 3 that consists of iterating between the Poisson and Vlasov equations and can be used even for cases where the Vlasov equation does not have an analytical solution. This approach will serve as the bedrock for the more elaborate cases that we have studied in this thesis, e.g., when the TEC comprises of grids (Chapter 4) or includes low-pressure plasmas (Chapter 5). As a way of efficiently obtaining a numerical solution, we present a particle tracing approach to calculate the electron density. Although we present the 1-dimensional case, this strategy can be expanded to include higher dimensions where an analytical solution to the Vlasov equation does not exist. This is an important issue, since in the absence of an analytical solution, the Vlasov equation lives in the 6-dimensional phase space, and thus finding a numerical solution to it can be extremely challenging and time-consuming. Our approach can also incorporate the analytical solution of the Vlasov equation; the 26 applicability of the model in this case is significantly improved by employing an asymptotic expansion to calculate the ill-behaving functions involved in the solution of the Vlasov equation. Dugan et al. [73] employed Langmuir’s theory of space charge in the presence of back-emission using a method similar to that of Hatsopoulos, i.e., starting with a current density and calculating the voltage that corresponds to this current. They developed an iterative approach, since the total current density depends on the unknown in the problem, which is the applied voltage. This approach naturally faces some of the challenges associated with Hatsopoulos' method, described above. It is indispensable to reliably include the presence of back-emission in the analysis of the emergent thermionic devices; by employing nanotechnology, ever smaller inter-electrode distances are achievable, which could lead to a higher temperature of the collector. In Chapter 3, after describing the model, we develop the physics of the TEC in the presence of back-emission and present the output current-voltage characteristics. The model presented in Chapter 3 can be used to evaluate the output characteristics of the device for a wide range of parameters, unless quantum tunneling (occurring at nanoscale inter-electrode distances) or relativistic effects (occurring at exceedingly high fields that are not relevant in thermionic converters) are prevalent. The next step is to model the behaviour of TECs in the presence of a positively-biased grid. As explained in the beginning of this section, it is established that the presence of grids can ameliorate the space charge effects. However, the model that has so far been used to analyze this situation is limited and may lead to spurious results. In Chapter 4, the physics of TECs in the presence of a grid is developed. It is demonstrated that a unique phenomenon of momentum gaps could arise in the phase space of electrons as they propagate in the inter-electrode distance. When one hears of forbidden values of electron momentum, one immediately thinks of the quantum mechanical effect arising from boundary conditions and the reflection of electron waves, and not a classical effect. Surprisingly, the momentum gap we report in Chapter 4 has a classical origin. The parallel is scientifically intriguing: while the quantum mechanical momentum gap can be viewed as being the result of having both forward-moving and backward-moving waves (leading to standing waves with quantized momenta) due to reflection from the boundary of a potential distribution, the classical momentum gap is a result of having both forward-moving and backward-moving particles due to reflection from the boundary of a potential well formed by space charge. We demonstrate that the highly nonlinear nature of the Vlasov-Poisson system leads to a subtle, yet critical, interplay between electron density and velocity distribution. The momentum gap affects this interplay in a fundamental manner: not accounting for this gap could prevent self-consistency in the Vlasov-Poisson system of equations. By noticing this momentum gap and including it in the transport problem, we are treating this interplay correctly for the first time. 27 The introduction of ions in the inter-electrode region has been the most ubiquitous method of overcoming the main challenges that have prohibited TECs to achieve widespread deployment. The positively charged ions neutralize the negative space charge of the electrons and also reduce the workfunction of the emitter by depositing the electropositive and neutral form of the ions on the surface of the emitter. However, in order for sufficient ion coverage on the surface, a high concentration of neutral atoms needs to be present due to the high temperature of the emitter. This is the basis of the operation of high pressure vapor diodes. However, the presence of a high concentration of neutral atoms in the inter-electrode region can cause scattering of the electrons and result in significant loss of electron current [14]. These high pressure diodes have been studied using phenomenological models or magnetohydrodynamics equations [13] [57]. Using nanotechnology, it is possible to engineer materials with inherently low workfunctions or to use a method such as intercalation to prepare much more stable surfaces with reduced workfunctions [74] [75] [76]. This allows for the operation of a thermionic converter in lower concentrations of neutral atoms. Therefore, it is possible to benefit from the advantages gained from the presence of ions without having a loss of current density due to collision of the electrons with neutral particles. In low-pressure plasmas, the mean free path (mfp) of the collision between the electrons and the neutral atoms is less than the inter-electrode distance, such that the presence of the neutral atoms may be neglected as a first approximation. For typical operating conditions of a TEC (e.g., a potassium reservoir temperature of 700 𝐾, leading to a vapor pressure of ~103 𝑃𝑎 [77]), the mfp between the charged particles (ions and electrons) and the neutral atoms is ~ 100 𝜇𝑚. The electron density is typically less than this in TECs and therefore a low-pressure plasma model is required for the analysis in this wide regime if the inter-electrode distance (𝑑) is on the order of ~ 100 𝜇𝑚 for the aforementioned temperature. The phenomenological model of the analysis of high pressure TECs relies on solving drift and diffusion equations, which in turn depend on the random walk motion induced by the collision with neutral particles. On the other hand, magnetohydrodynamic equations need to be truncated at some level (since the moments of the Boltzmann equation are solved for, rather than the transport equation itself [54]) and the diffusive (random walk) features of the neutral atoms’ movements are introduced. Therefore, different sets of equations apply to the operation of TECs in low and high pressure plasma cases. There are some works on the physics of low-pressure plasmas applied to thermionic converters from several decades ago, notably due to Auer [78] and Hatsopoulos [57]. Both of these works are based on finding analytical solutions to the Vlasov-Poisson system and have a limited range of applicability and predictive power [3]. 28 In Chapter 5, we adopt a different strategy based on a direct solution of the Boltzmann transport equation in the form of a Vlasov-Poisson system of equations, but instead of trying to find an analytical solution for the system, we adopt a numerical, iterative approach between the Poisson and Vlasov equations. This allows us to solve the Vlasov equation even in cases that an analytical solution to the Vlasov-Poisson system does not exist or becomes intractable (in the presence of several peaks in the motive). We include all the electron-electron, electron-ion and ion-ion interactions in the Poisson equation in a mean-field fashion and introduce the resulting motive as a force term inside the Vlasov equation. We further demonstrate that under certain initial conditions, momentum gaps arise in the phase space of the electrons and ions, similar to the case where a grid is present in the inter-electrode region of a TEC (in the absence of ions). Some recent work on the analysis of high pressure plasmas in the presence of high electric fields have been done by Go et al. [79] [80] . The fields are high enough to induce breakdown of the gas in the inter-electrode region [81]. The analysis that we provide in this thesis is different in the sense that we are interested in the low-pressure regime of operation. We will illustrate that, depending on the flux of the ions, different scenarios could arise. At sufficiently low ion fluxes, the motive barrier is lowered, but the overall shape of the motive remains similar to the case where no ions are present. However, at higher ion fluxes, plasma oscillations could arise and the motive is significantly altered. We finally extend the model to include electron and ion emission from both electrodes. The experimental results of this thesis are presented in the appendices. In Appendix A, the model developed in Chapter 2 is used to solve the reverse problem starting from the experimental I-V characteristics that are obtained from a MWCNT-based LITEC. The parameters of the device are then extracted from the I-V characteristics. In Appendix B, the reduction of the workfunction of the potassium intercalated MWCNT in an in-situ environment under LITEC conditions is experimentally studied and the results of Chapters 2 and 5 are used to extract the amount by which the workfunction of the MWCNTs is reduced. It is crucial to study the reduction of the workfunction of CNTs under in-situ conditions, since this behaviour is highly temperature dependent and, therefore, an ex-situ method of workfunction measurement such as ultraviolet photoemission spectroscopy would not be sufficient. 29 Chapter 2- An analytical approach to the analysis of thermionic energy converters2 The existing models, which were developed for macroscopic converters, are not capable of describing all aspects of nanostructured devices. In this chapter, we present extensions on the current method to evaluate the output characteristics of thermionic converters with a higher precision than the existing models and the ability to simulate a broader range of parameters, including temperatures, active surface areas, inter-electrode distances, and workfunctions. These features are crucial for the characterization of emergent devices due to the unknowns involved in their internal parameters; the model's high numerical precision and flexibility allows one to solve the reverse problem and to evaluate the internal parameters of the device from a set of simple experimental data. Theory and model It is assumed that electrons are thermionically emitted from the cathode as governed by the Richardson-Dushman equation (equation (1-1)), 𝐽𝑠𝑎𝑡 = 𝐴 𝑇𝐸2 exp(−𝜙𝐸/𝑘𝐵 𝑇𝐸), where 𝐽𝑠𝑎𝑡 denotes the saturation current density, i.e. the maximum current density attainable from the emitter at the emitter’s temperature, 𝑇𝐸, 𝑘𝐵 is the Boltzmann constant, and 𝐴 = 1.202 × 106 𝐴𝑚−2𝐾−2 is the Richardson-Dushman constant [82]. The electrons reaching the other electrode (collector or anode) are completely collected and run through the external circuit. The methodology in this section follows the earlier description of operation of TECs in the space charge limited regime [57], [29]. It is described here how the numerical accuracy and the range of applicability can be further enhanced using a new algorithm. Moreover, the I-V characteristics of the TEC are obtained for different operating regimes of device operation. The high numerical accuracy, applicability to a wide range of parameters, and modelling of the entire I-V characteristics, in addition to a stable algorithm, make solving the more challenging reverse problem possible. The general method to solve the direct problem is to couple Poisson’s and Vlasov’s equations. Vlasov’s equation applies when electrons in the inter-electrode region are approximated as a collisionless ideal gas [83]. In the simplest one-dimensional case, i.e., when a uniform electric field is applied only in the direction normal to the emitter 2 This chapter is reproduced from [2] with minor modifications. 30 and collector’s surface (e.g. along the x axis), an analytical solution is possible. The solution includes functions and integrals that can be solved numerically. The operation of the TECs in the space charge mode is discussed in Ref. [57]. We provide a summary of the key points along with important equations for the sake of clarity. We also propose an algorithm that can be utilized to increase the accuracy of the method proposed in Ref. [14] to calculate the I-V characteristics. Figure 2-1. (a) A typical I-V curve with the corresponding motive diagrams in each operation regime. (b)-(c) Typical motive diagrams corresponding to different operation regimes of a vacuum thermionic converter under negative (b) and positive biases (c). The subscripts “E” and “C” denote emitter and collector, respectively. 𝜓, 𝜇 and 𝜙 denote the motive, electrochemical potentials and workfunctions, respectively. 31 A TEC can operate in two or three distinct regimes (or modes), depending on the electric potential difference between the electrodes (Figure 2-1). The width of each mode and the onset of the subsequent mode depend on the internal parameters of the TEC, such as temperature, inter-electrode gap distance, and the workfunctions. Figure 2-1.(a) shows the I-V characteristics of a typical TEC in its various operating regimes. The appropriate motive diagrams that are the analogs of the band-diagrams in solid state physics are also included. The emitter and collector are assumed to have the same Richardson-Dushman constant and maintained at 𝑇𝐸 and 𝑇𝐶 , respectively. The motive diagram of a TEC operating under negative and positive biases is depicted in Figure 2-1. (b-c). At sufficiently negative biases, where the chemical potential of the emitter lies considerably below that of the collector, the device is in the retarding mode. The electric field opposes the transport of electrons from the emitter toward the collector and only a fraction of them with sufficient kinetic energy to overcome the electrostatic barrier can make it to the other side. Therefore, the current in this regime can be calculated as 𝐽 = 𝐴 𝑇𝐸2 exp(−(𝜙𝐶 + 𝑒𝑉)/𝑘𝑇𝐸). It is noted that in this mode, due to the low concentration of electrons, the electron-electron interaction is negligible and, therefore, not included in the calculations. In other words, all electrons are assumed to be subject to the same opposing electric field in the inter-electrode gap distance; the spatial changes of the field due to the presence of electrons in the inter-electrode region are not considered. As the concentration of electrons rises by lowering the potential energy of the collector with respect to that of the emitter, the electric field arising from the presence of electrons in the inter-electrode region must be taken into account. The electrons ejected from the emitter effectively face a higher local potential barrier due to the electrostatic repulsion arising from the other electrons. This characterizes the space charge limited regime. As observed from Figure 2-1, the presence of the electrons in this region results in a peak in the inter-electrode potential (𝜓 = 𝜓𝑀 at position 𝑥 = 𝑥𝑀). Electrons situated between the surface of the emitter and 𝑥𝑀 can have both positive and negative velocities. Only those electrons with kinetic energy equal to or greater than 𝜓𝑀 (motive of the point just outside the emitter is set as the origin) are able to pass over the maximum potential. On the other hand, once electrons overcome 𝜓𝑀, their potential energy only decreases for values of 𝑥 higher than 𝑥𝑀. Consequently, they can only drift toward the collector. As the applied bias increases, 𝑥𝑀 moves consistently toward the emitter and, eventually, 𝜓𝑀 does not take place in the inter-electrode gap distance. The device has entered the saturation mode. The value of current density in this case is equal to 𝐽𝑠𝑎𝑡 and remains constant by increasing the electric field. This approximation is valid for relatively low applied biases, where the Schottkey effect becomes important at 32 higher biases, where the effective workfunction is lowered according to 𝜙𝑒𝑓𝑓 = 𝜙𝐸 − √𝑒3 𝐹4 𝜋𝜖0, where 𝜖0 is the vacuum constant and 𝐹 is the magnitude of the electric field [84][85]. Space charge limited regime The space charge limited regime is characterized by two boundary points. The onset of the space charge limited regime, where 𝜓𝑀 takes place at the point just outside the collector, is known as the critical voltage. This point is distinguished by 𝐽𝐶 and 𝑉𝐶 as the critical current density and applied voltage, respectively. On the other hand, when 𝜓𝑀 coincides with the point just outside the emitter, the device is at the saturation point. This point is identified by 𝐽𝑠𝑎𝑡 and 𝑉𝑠𝑎𝑡. In order to obtain the overall I-V characteristics, the saturation current and voltage are determined. Then, as explained in the following section, the critical point is evaluated. The values of the current and voltage between these two boundaries are subsequently computed from coupling Vlasov’s and Poisson’s equations. Starting with a current, Vlasov’s equation outputs the velocity distribution function, 𝑓(𝑥, 𝑣𝑥). The electron density as a function of 𝑥, 𝑛(𝑥) is evaluated by integrating 𝑓(𝑥, 𝑣𝑥). The resulting 𝑛(𝑥) is plugged into Poisson’s equation. The solution to Poisson’s equation can be evaluated numerically to obtain the bias particular to that current. Theory of the space charge limited regime In order to calculate the characteristics in the space charge limited mode, the electron velocity distribution function, 𝑓(𝑥, 𝑣𝑥), the electron number density, 𝑛(𝑥), and the dimensionless Poisson’s equations are obtained and solved numerically. In the case of collisionless electrons, Vlasov’s equation can be solved to determine the distribution function of electrons. At 𝑥𝑀, the distribution function, 𝑓(𝑥𝑀 , 𝑣𝑥) takes the form of a half-Maxwellian at temperature 𝑇𝐸; i.e., for 𝑣𝑥 > 0, the distribution is Maxwellian, whereas for 𝑣𝑥 < 0, the distribution is zero. 𝑣𝑥 represents the velocity of the electron along the x direction. The equation for 𝑓(𝑥𝑀 , 𝑣𝑥) can be used as the boundary condition to calculate the distribution function. Within the region 𝑥 > 𝑥𝑀, in Figure 2-1, the electrons are accelerated by the negative space charge since the gradient of the motive is always negative. The minimum value of the velocity at point 𝑥, i.e., 𝑣𝑥,𝑚𝑖𝑛 is therefore equal to the velocity component of electrons starting with zero initial velocity at 𝑥𝑀. Consequently these electrons undergo a motive difference of 𝜓𝑀 − 𝜓(𝑥) corresponding to 𝑣𝑥,𝑚𝑖𝑛 = (2 𝜓𝑀−𝜓(𝑥)𝑚)1/2and therefore, when 𝑥 > 𝑥𝑀 the solution to the Vlasov equation is equal to 33 𝑓(𝑥, 𝑣𝑥) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘 𝑇𝐸)32𝑒𝑥𝑝 (𝜓𝑀−𝜓(𝑥)𝑘𝑇𝐸−𝑚𝑣𝑥22 𝑘𝑇𝐸) ×𝛩(𝑣𝑒𝑥 − 𝑣𝑥,𝑚𝑖𝑛), (2-1a) where 𝑚 is the mass of the electron, 𝑣 = (𝑣𝑥2 + 𝑣𝑦2 + 𝑣𝑧2)1/2 is the electron speed, and Θ is the Heaviside step function. On the other hand, at 𝑥 < 𝑥𝑀, the electrons are decelerated by the electric field and can have both positive and negative velocities along the x direction. This leads to a negative minimum velocity along the x direction, 𝑣𝑥,𝑚𝑖𝑛 = −(2 𝜓𝑀−𝜓(𝑥)𝑚)1/2. Therefore, when 𝑥 ≤ 𝑥𝑀, the velocity distribution function is given by 𝑓(𝑥, 𝑣𝑥) = 2𝑛(𝑥𝑀) (𝑚2 𝜋𝑘 𝑇𝐸)32exp(𝜓𝑀 − 𝜓(𝑥)𝑘𝑇𝐸−𝑚𝑣𝑥22 𝑘𝑇𝐸)× Θ(𝑣𝑥 + 𝑣𝑥,𝑚𝑖𝑛). (2-1b) Based on this explanation, electrons on either side of 𝑥𝑀 also undertake a half-Maxwellian velocity distribution with the peak shifted from 0 to 𝑣𝑥,𝑚𝑖𝑛. The next step is to calculate the electron number density, 𝑛(𝑥), by integration of the velocity distribution function over the entire velocity space at 𝑥. This results in 𝑛(𝑥) = 𝑛(𝑥𝑀) exp(𝛾) {1 − erf(𝛾),1 + erf(𝛾) ,𝑥 > 𝑥𝑀𝑥 ≤ 𝑥𝑀, (2-2) where 𝛾 ≡ (𝜓𝑀 − 𝜓(𝑥))/𝑘𝐵𝑇𝐸 is the dimensionless potential and 𝑒𝑟𝑓 is the error function. In order to obtain the dimensionless Poisson’s equation, the variable 𝑥 in the Poisson equation (1-6) is divided by the Debye length given by 𝐿D ≡ (𝜖0𝑘𝐵𝑇𝐸2 𝑒2𝑛𝑒(𝑥𝑀))12 [83]. The resulting equation is 𝑑2𝛾𝑑𝜁2= exp(𝛾) × {1 + erf(𝛾1/2),1 − erf(𝛾1/2) ,𝜁 ≥ 0𝜁 < 0, (2-3) 34 where 𝜁 ≡ (𝑥 − 𝑥𝑀)/𝐿D is the dimensionless distance. Noting that 𝜓𝑀 is a stationary point of 𝜓(𝑥), the corresponding boundary conditions are obtained as 𝛾(𝜁=0) = 0 and (𝑑𝛾𝑑𝜁)(𝜁=0)= 0. Double-integrating the equation and utilizing the boundary conditions leads to 𝜁 = − ∫𝑑𝛼[exp(𝛼) + exp(𝛼) erf (𝛼12) − 2 (𝛼𝜋)12− 1]1/2𝛾0 (2-4a) for 𝜁 < 0 and 𝜁 = ∫𝑑𝛼[exp(𝛼) − exp(𝛼) erf (𝛼12) + 2 (𝛼𝜋)12− 1]1/2𝛾0 (2-4b) for 𝜁 > 0. It is noted that here the term Debye length is used in accordance with the definition in plasmas and electrolytes and is not related to the other contributions that the Dutch physicist Peter Debye has made to other areas of physics (e.g., the Debye model of phonons). This length also appears in Thomas-Fermi screening of the electric field due to the interactions between electrons and nuclei in a solid [59]. The integral on the right hand side of equation (2-4) was solved numerically with higher precision and with 3 orders of magnitude of wider range compared to [14], as plotted in Figure 2-2 (a). As depicted in Figure 2-2 (b), having a limited range for the dimensionless variables and using interpolation as prescribed in ref. [14], can lead to substantial errors. 35 Figure 2-2.(a) The dimensionless motive (𝛾) and dimension distance (𝜁) plotted as the solution of the coupled Poisson and Vlasov equations. (b) The values of the dimensionless motive, 𝛾 are fitted using linear, quadratic and cubic polynomial in the range 0 < 𝜁 < 15, as tabulated in Ref. [14]. These fitted polymonials are used to extrapolate the values of 𝛾 for 𝜁 > 15. As observed, having a limited range for dimensionless motive leads to substantial errors in case that a limited range of the values of the dimensionless variables are used. A high precision and a wide range for the values of the dimensionless motive and distance are imperative in order to be able to solve the reverse problem and extract the device parameters from the experimental I-V curves. The calculation of the values of the integrals in equation (2-4) requires dealing with functions that change drastically depending on the values of 𝛼. Furthermore, the limits of integration involve large numbers (the values of 𝛾), growing ever larger once used as inputs to exponential functions. The values of 𝜁 for 0 < 𝛾 < 14 are tabulated in Ref. [14]. The existing strategy was to extrapolate these values to the desired specifications of the device. However, as will be shown in the following section, the extrapolated values obtained by fitting to different polynomials do not lead to desirable accuracy at higher values of 𝜁. The greatest error arises when determining the values of the critical current and voltage. In order to evaluate the integrals, the variable-precision arithmetic feature from MATLAB’s Symbolic Math Toolbox was employed to compute each element to at least 7 significant digits of accuracy. After calculating the integrals numerically, it was noticed that for 𝜁𝐶 > 10, 𝛾(in the range 𝛾 > 10.94) can be fitted to a 3-term function in the form 𝛾𝐶 = 𝑎 𝜁𝐶𝑏 + 𝑐 + 𝑒 𝜁𝐶𝑑, where 𝑎 = 0.7388, 𝑏 = 1.333, 𝑐 = 2.198, 𝑑 =0.6712 and 𝑒 = −1.476, with a root mean square error (RMSE) of 1.144. Moreover, for 𝜁𝐸 > 2, 𝜁𝐸 = 𝑎 exp(−𝑏𝛾𝐸) + 𝑐, where 𝑎 = −1.497, 𝑏 = 0.5125, 𝑐 = 2.554, and 𝑅𝑀𝑆𝐸 = 0.0004346. 36 Figure 2-3. Comparison of the estimated value of the critical point voltage, 𝑉𝐶, when calculated using the assumptions made by Hatsopoulos, 𝑉𝐶,𝐻𝑎𝑡, and the actual value of 𝑉𝐶 obtained from iteratively solving the problem at different (a) inter-electrode distances and (b) emitter temperatures. These results are calculated 37 for typical values of device parameters of 𝜙𝐸 = 4.6 𝑒𝑉, 𝜙𝐶 = 4.0 𝑒𝑉. The radius of the heat spot, 𝑟, is assumed to be 100 𝜇𝑚. In part (a), 𝑇𝐸 = 2000 𝐾 and in part (b), 𝑑 = 700 𝜇𝑚. The value of the current density can be calculated by integrating 𝑒𝑓(𝑥, 𝑣𝑒)𝑣𝑒𝑥 at any arbitrary 0 < 𝑥 <𝑑, (𝑑 being the inter-electrode distance) over the entire velocity space. The result is 𝐽 =𝑒 𝑛(𝑥𝑀) (2𝑘𝐵𝑇𝐸𝜋𝑚)1/2. This equation can be used to replace 𝑛(𝑥𝑀) in the definition of Debye length and directly relate it to the current as 𝐿𝐷 = (𝜖02𝑘𝐵32 𝜋𝑚 𝑒2)14 𝑇𝐸3/4𝐽1/2. The final step is to link 𝐽 to the maximum current density obtainable from the emitter, i.e., 𝐽𝑠𝑎𝑡 by calculating the number of electrons with 𝑣𝑥>0 at the point just outside the emitter. In other words, 𝐽 = 𝐽𝑠𝑎𝑡 when 𝑥𝑀 = 0. This leads to 𝐽𝑠𝑎𝑡 = 𝐽 exp(𝛾𝐸), where 𝛾𝐸 = 𝛾(𝑥 = 0). Based on these equations and values from Figure 2-2.(a), the following procedure is proposed to obtain the I-V characteristics of the entire device. Current-voltage characteristics at the critical and saturation points At the saturation point, the maximum motive takes place at the point just outside the emitter (𝑥𝑀 = 0 and 𝐽 = 𝐽𝑆). Therefore, 𝜁𝐸 = 0 and 𝜁(𝑥 = 𝑑) = 𝜁𝐶 =𝑑−0𝐿𝐷=(2 𝜋𝑚 𝑒2𝜖02𝑘𝐵3 )14 𝐽𝐸𝑠1/2𝑑𝑇𝐸3/4 . Based on this value of 𝜁𝐶, 𝛾𝐶 is interpolated from Figure 2-2(a). Following the notations in Figure 2-3, with the origin of the x coordinate at the point just outside the emitter, the saturation voltage, 𝑉𝑠𝑎𝑡 = 𝑉𝑎𝑝𝑝, is calculated as (𝜙𝐸 − 𝜙𝐶 −𝛾𝐶𝑘𝑇𝐸). In calculating the critical point characteristics, it is noted that at this point, the maximum motive occurs just outside the collector (𝑥𝑀 = 𝑑). Therefore, 𝛾𝐸 = ln(𝐽𝑠𝑎𝑡/𝐽𝐶), |𝜁𝐸| = (0−𝑑𝐿𝐷) = (2 𝜋𝑚 𝑒2𝜖02𝑘𝐵3 )14 𝐽𝐶1/2𝑑𝑇𝐸3/4 , and 𝜁𝐶 = 0. In order to compute 𝐽𝐶, Hatsopoulus made the assumption that at this point, the distance |𝜁𝐸| reaches its maximum value (2.554 from Figure 2-2(a)). Therefore, 𝐽𝐶 can be calculated from this value of |𝜁𝐸| and then replaced into the Richardson-Dushman equation to evaluate the critical point voltage, 𝑉𝐶. This estimation, as presented in Figure 2-3 (a-b), leads to numerical errors. The value of the 𝑉𝐶 calculated using the assumption |𝜁𝐸|=2.554 is denoted by 𝑉𝐶,𝐻𝑎𝑡. These values are compared with the 𝑉𝐶 values that are obtained using the following proposed method. According to Figure 2-3, the approximation made can lead to considerable error in determining the value of 𝑉𝐶. We propose this method to calculate the exact value of 𝑉𝐶 as follows: 38 a. As a first approximation, 𝐽𝐶 and 𝑉𝐶 are calculated following the method explained above. b. Values of 𝛾𝐸 = ln(𝐽/𝐽𝑠𝑎𝑡) are computed for 10−6 <𝐽𝐽𝐶< 103. c. 𝜁𝐸 is computed for each 𝛾𝐸, based on their relation from Figure 2-2.(a). d. 𝜁𝐶 is calculated from 𝜁𝐶 =𝑑−𝑥𝑀𝑥0=𝑑𝑥𝑜− 𝜁𝐸. If 𝜁𝐶 < 0, the maximum motive occurs outside the inter-electrode gap and, therefore, the device does not enter the space charge regime, and operates in the retarding mode. The smallest current that sets 𝜁𝐶 = 0, is equal to the actual value of 𝐽𝐶. e. 𝑉𝐶 is calculated from the Richardson-Dushman equation in the retarding regime: 𝐽 =𝐴 𝑇𝐸2 exp(−(𝜙𝐶 + 𝑒𝑉)/𝑘𝑇𝐸) The value of 10−6 for the ratio of currents allows for calculations of actual critical point currents up to 6 orders of magnitude different from the estimated value. If the integrals are not calculated for a wider range, it will not be possible to calculate the exact value of 𝐽𝐶, as large values of 𝛾𝐸 are required in Step c. Output current-voltage characteristics in the space charge limited regime In order to calculate the I-V characteristics, the following algorithm is proposed: a. A procedure similar to the Steps b-d in calculating the critical point is employed, except that the current is chosen within the range of 𝐽𝐶 < 𝐽 < 𝐽𝑠𝑎𝑡. b. The values of 𝛾𝐶 are computed from the relation between 𝛾𝐶 and 𝜁𝐶, as illustrated in Figure 2-2(a). c. The output voltage is calculated from Figure 2-1 as (𝜙𝐸 − 𝜙𝐶) + 𝑘𝐵 𝑇𝐸(𝛾𝐸 − 𝛾𝐶). 39 Figure 2-4. The entire I-V characteristics and the important boundary points (critical and saturation voltages and currents) were calculated for a wide range of device parameters. The effects of altering (a) the inter-electrode gap size and (b) the hot spot temperature on the I-V characteristics are included in here. The values of the currents and voltages at the boundary points (critical and saturation points, represented by subscripts “C” and “sat”, respectively), are tabulated in here. For both (a) and (b), 𝜙𝐸 = 4.6 𝑒𝑉 and 𝜙𝐶 = 4.0 𝑒𝑉, 𝑟 = 100 𝜇𝑚. For part (a), 𝑇𝐸 = 2,000 𝐾 and for part (b), 𝑑 = 500 𝜇𝑚. 40 Figure 2-4.(a) and Figure 2-4.(b) demonstrate the effects of changing the gap distance and the emitter’s temperature on the I-V characteristics of the device, respectively. In the retarding and saturation regimes, the current is determined using the Richardson-Dushman equation with the proper value of the energy barrier that the electron encounters when travelling from the emitter to the collector. On the other hand, in the space charge limited regime, numerical solutions of equation (2-2) are used to calculate the characteristics. The entire I-V characteristics are eventually obtained by combining the I-V curves within these three regions. Non-uniform temperature profile at the emitter So far, we have assumed that the temperature is uniform across the illuminated area. This can be a reasonable assumption in the case of the bulk TEC device with large dimensions. However, a detailed analysis of the new generation of TECs with much smaller dimension reveals the presence of temperature gradients. For instance, in a real LITEC device, inspection of the temperature profile reveals an average temperature gradient of several degrees of kelvin over a micrometer at a typical input laser power of 50 mW [86]. Therefore, a uniform temperature profile is not realistic and the influence of the existence of different areas with various temperatures has to be included in the model. In order to model this effect, the circular heat-spot area is divided into different rings with different temperatures (𝑇(𝑟𝑛)), with 𝑟𝑛 being the distance of the nth ring from the center of the beam. Each ring is considered as a thermionic converter contributing to the overall current. Linear and Gaussian temperature profiles were used to investigate this effect. In the case of a linear temperature profile, 𝑇(𝑟𝑛) = (𝑇0 − 𝑇𝑚𝑎𝑥)𝑟𝑟0+ 𝑇𝑚𝑎𝑥, this leads to an average temperature of 𝑇𝑎𝑣𝑒𝑟 =1𝜋𝑟02 ∫ 𝑇(𝑟𝑛) 2 𝜋𝑟 𝑑𝑟 =13(2 𝑇0 + 𝑇𝑚𝑎𝑥𝑟00), where 𝑇𝑚𝑎𝑥 is the maximum temperature, 𝑇0 = 𝑇(𝑟0), and 𝑟0 is the radius of the thermionic emission area (Figure 2-5.(a)). The value of 𝑇0 is chosen such that, the contribution from the ring at 𝑇 = 𝑇0, is less than 10−4 % to the overall current. In the case of a Gaussian temperature distribution, 𝑇(𝑟) = 𝑇𝑚𝑎𝑥 exp (𝑟2𝛼𝑟02), where 𝛼 = ln(𝑇0/𝑇𝑚𝑎𝑥), leads to 𝑇𝑎𝑣𝑒𝑟 = 𝑇𝑚𝑎𝑥𝛼(−1 + exp(1/𝛼)). The contributions of different temperatures when the thermionic emission area is divided into 10 rings are represented in Figure 2-5.(b). As observed, the average spatial temperature of the spot (2,020 K) does not fully represent the thermionic emission behaviour of the spot. The reason is the exponential dependence of the emission current on the temperature of the emitter. Therefore, the defining temperature of the device closer to the maximum temperature than the spatial average temperature. Remarkably, one single temperature cannot fully represent the system in any case. At negative biases, higher-temperature elements with less area contribute predominantly to the current (2,710 𝐾 with an area of 5.59 × 103 𝜇𝑚2). On the other hand, near the saturation regime, the higher-41 surface-area rings with relatively lower temperatures dominate. The reason is the high- and low-temperature areas facing the same energy barrier at sufficiently negative biases (𝜙𝐶 + 𝑒|𝑉|). However, in the space charge limited regime, higher-temperature sectors feel a larger barrier due to the increased space charge, since they have a higher current density. 42 Figure 2-5. (a) Temperature profile (assuming linear temperature distribution) on the hot spot of the CNT-based LIETC. (b) Contributions of different temperature contour rings to the total current across a hot spot with a radius, 𝑟 = 200 𝜇𝑚. Summary A method was proposed based on the simultaneous solutions of Vlasov and Poisson equations to calculate the output characteristics of thermionic converters. The numerical accuracy and the range of applicability of the method were substantially enhanced by applying robust algorithms. It was demonstrated that this method can be employed to solve the reverse problem and calculate the internal parameters, e.g. workfunctions, inter-electrode distances, effective areas and temperatures of the thermionic converters. It was shown that the emergent nanomaterial-based thermionic, though operating on the same overall principles, have some distinct features such as spatially varying workfunctions and the presence of temperature gradients, which require more sophisticated simulation methods to study them. Due to the high stability of the proposed method, a broad range of parameters can be solved for and, therefore, the method can address the numerical difficulties associated with new thermionic devices. As an example of its application, this method was used to extract the parameters of a light induced thermionic conversion device based on carbon nanotubes (Appendix A). 43 Chapter 3- A self-consistent approach to the analysis of thermionic devices3 As detailed in the previous Chapters, thermionic energy converters are commonly modelled using the Vlasov-Poisson system of equations under various limitations and approximations. With the ever-growing demands of emergent thermionic devices, more comprehensive approaches are needed in order to be able to treat a broader range of device configurations and operational parameters. Here, we propose a self-consistent approach that, by iterating between the Poisson and Vlasov equations, does not rely on the existence of an analytical solution to the latter. Specifically, we present a particle-tracing implementation of this method for solving the system numerically in an efficient manner. In the case where an analytical solution does exist, we present an asymptotic expansion of the ill-behaving functions that arise; this approach improves the effectiveness of the method in the deep space charge mode. We also demonstrate the applicability of this approach in the presence of back-emission. Theory and model In the first section, we develop the model in the absence of back-emission and show that the results are in agreement with our improved-and-extended-Langmuir (IEL) solution reported earlier [2] (Chapter 2). In the subsequent sections, we use this model to analyze a TEC with substantial back-emission from its collector. It is assumed that the electron flux from the electrodes into the inter-electrode space is determined by Richardson-Dushman equation. However, if this equation is not applicable under some circumstances (for instance from lower dimensional nanomaterials such as graphene [87]), the proposed model in this chapter is still equally valid, since in principle, the Richardson-Dushman equation serves as a boundary condition for the particle flux, and can be replaced by the appropriate equation that describes the material. 3 This chapter is reproduced from [3] with minor modifications. 44 Self-consistent solution in the absence of back-emission In a 1-D case and in the absence of magnetic fields, the steady-state collisionless Vlasov equation (1-3) is used [18] [54]. It is assumed that the electrons adopt a hemi-Maxwellian distribution (equation 1-5) at the point of maximum motives and follow the statistics governed by equations (2-1a) and (2-1b) at the points lying on the right and left hand side of the maximum motive, respectively. In the next section, we will provide arguments on the validity of this assumption in the case of the emergent TECs and also in the presence of back-emission. The maximum motive, 𝜓𝑀, coincides with the points just outside the emitter or collector at the start of the flowchart in Figure 3-1, in positive and negative applied voltages, respectively. However, the maximum motive can move inside the inter-electrode region in subsequent iterations. In that case, the device is operating in the SCL mode. If the maximum motive coincides with the points just outside the emitter or collector at the end of the loop, the device is operating in the saturation or retarding modes, respectively. As we will show, our strategy can be used in these modes as well. 45 Figure 3-1. Flowchart representation of the proposed self-consistent approach to calculate the output characteristics of TECs. Initially, the electron density of the electrons is assumed to be zero, and therefore only the Laplace equation is solved. In later iterations, the Poisson and Vlasov equations, or alternatively, Poisson and the particle tracing equations are self-consistently solved to reach convergence. The velocity distribution function, 𝑓(𝑥, 𝑣), can be calculated by considering the ranges of applicable velocities along different points in the inter-electrode region as detailed in Chapter 2. This range of velocities is again due to the fact that electrons originate from the emitter, and the ones that possess sufficient kinetic energy to overcome the 𝜓𝑀 barrier, only move to the right at 𝑥 > 𝑥𝑀. The equations governing the velocity distribution function of a typical TEC in the absence of back emission and in the 1-D case were derived in Chapter 2 (equations (2-1)). The resulting electron density, 𝑛(𝑥), is calculated by integration of the velocity distribution function as 𝑛(𝑥) = ∫ ∫ ∫ 𝑓(𝑥, 𝑣) 𝑑𝑣𝑥 𝑑𝑣𝑦 𝑑𝑣𝑧∞−∞∞−∞∞−∞ as presented by equation (2-2). The electron density is subsequently substituted into Poisson’s equation to yield 46 𝑑2𝜓(𝑥)𝑑𝑥2= −𝑒2𝜖0𝑛(𝑥𝑀) exp(𝛾𝐸) {1 − erf(𝛾𝐸1/2),1 + erf(𝛾𝐸1/2 ) ,𝑥 > 𝑥𝑀𝑥 ≤ 𝑥𝑀, (3-1) subject to the boundary conditions, 𝜓(0) = 0 and 𝜓(𝑑) = 𝑒 𝑉𝑎𝑝, where 𝑒 is the electron’s charge (negative value), 𝜖0 is the permittivity of vacuum, and 𝑉 is the voltage difference between the collector and the emitter in the case that the workfunctions of the emitter and the collector are equal. In the more general case where the workfunctions are not equal, the electric potential difference between the points just outside the collector and the emitter is equal to 𝑉𝑎𝑝 − (𝜙𝐸 − 𝜙𝐶)/𝑒, where 𝑒𝑉𝑎𝑝 is the difference between the Fermi levels of the collector and the emitter (the applied voltage). Also, 𝛾𝐸 ≡ (𝜓𝑀 − 𝜓(𝑥))/𝑘𝐵𝑇𝐸 (cf. the dimensionless Poisson equation (2-3)). Equation (3-1) can be solved at each iteration by noting that 𝜓(𝑥) can be written as the sum of two functions, 𝜓𝑙(𝑥) + 𝜓𝑑𝑝(𝑥), where 𝜓𝑙(𝑥) represents the solution to the Laplace equation, 𝑑2𝜓𝑙(𝑥)𝑑𝑥2= 0, with 𝜓𝑙(0) = 0 and 𝜓𝑙(𝑑) = 𝑒 𝑉 boundary conditions; 𝜓𝑑𝑝(𝑥) is the solution to the Poisson equation with the homogeneous Dirichlet boundary conditions, namely, equation (3-1) subject to 𝜓𝑑𝑝(0) = 𝜓𝑑𝑝(𝑑) = 0, as the boundary conditions. In this 1-D example, 𝜓𝑙(𝑥) =𝑥𝑑𝑉. 𝜓𝑑𝑝(𝑥) is solved numerically by means of spatial discretization of the Laplacian operator, 𝑑2𝑑𝑥2≈1Δ𝑥2[ −2 1 0 0 01 −2 1 0 0000100⋱⋱0⋱−2101−2 ] , (3-2) and taking the inverse of the resulting matrix. Δ𝑥 represents the discretization distance. Since the sampling points remain the same within the entire calculations for a given set of variables, the inverse of the Laplacian matrix can be calculated only once to save on computation cost. Or more efficiently, it is possible to use sparse matrices to solve this equation. The sparse matrix functions of MATLAB were used to this end. Stephanos and Meir et al. [15] [14] arrived at an equivalent form of equation (3-1) and used a non-linear solver to calculate the motive in the inter-electrode region.) An important numerical difficulty in solving equation (3-1) can arise at higher values of 𝛾𝐸 at 𝑥 > 𝑥𝑀 where the error function approaches 1. This issue can be bypassed by using an asymptotic expansion of the error function [88] and substituting the resulting function in equation (3-1): 47 exp(𝛾𝐸) (1 − erf(𝛾𝐸1/2))=1√𝜋(𝛾𝐸−1/2 −12𝛾𝐸−32 +34𝛾𝐸−52 − ⋯+(−1)𝑛+11 × 3 × … (2𝑛 − 3)𝛾𝐸𝑛−1/22𝑛−1+ ⋯ ). (3-3) Using the first 6 terms in equation (3-3) leads to numerical errors less than 10−4 % for values of 𝛾𝐸 >20. After Poisson’s equation is solved, the electric potential is updated as a mixture of the previous solution and the new solution to avoid large jumps in the potential, i.e., 𝜓𝑑𝑝,𝑛(𝑥) = (1 − 𝛿)𝜓𝑑𝑝,𝑛−1(𝑥) +𝛿𝜓𝑑𝑝(𝑥) , where 𝛿 is the mixing ratio, and the subscript 𝑛 represents the iteration number. This step is necessary due to the non-linearity of equation (3-1), leading to a strong dependence of electron density on the potential profile. The mixing ratio needs to be chosen small enough so that large oscillations do not occur in smaller values of 𝑛, where the solutions of the Poisson and Vlasov equations are strongly decoupled. Values of 𝛼 ≤ 0.1 were found to be suitable for the problems studied in this Chapter. The current density can be equal to the maximum saturation current density, 𝐽𝑠𝑎𝑡, determined by the Richardson-Dushman equation, if electrons do not face a potential barrier in the inter-electrode distance, i.e., when 𝑥𝑀 = 0. Therefore, 𝐽 = 𝐽𝑠𝑎𝑡 exp(−𝛾𝐸(𝑥 = 0)) = 𝐽𝑠𝑎𝑡 exp (−𝜓𝑀𝑘𝐵 𝑇𝐸), which can be used to calculate the current at each iteration. The current density at each iteration is equal to the integral of the product of the x component of the velocity and the velocity distribution function: 𝐽 = 𝑒 ∫ ∫ ∫ 𝑣𝑥𝑓(𝑥, 𝑣) 𝑑𝑣𝑥 𝑑𝑣𝑦 𝑑𝑣𝑧∞−∞∞−∞∞−∞= 𝑒 𝑛(𝑥𝑀) (2 𝑘𝐵 𝑇𝐸𝜋𝑚)1/2. (3-4) The result in equation (3-4) was calculated using equation (2-1). This value of current density is expectedly independent of position. The electron density at the point of maximum motive, 𝑛(𝑥𝑀), can be derived from this result. 48 The new motive is subsequently incorporated into the Vlasov equation and the loop is repeated until the motive converges. The convergence criterion is set such that the cumulative root mean square change in motive, √∑(𝜓𝑑𝑝(𝑥𝑖)−𝜓𝑑𝑝,𝑛−1(𝑥𝑖))2∑(𝜓𝑑𝑝,𝑛−1(𝑥𝑖))2 , is less than 1 %. Figure 3-2 presents the evolution of the motive as a function of the iteration number in the simulation. (The parameters of the TEC under study are 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 2,000 𝐾 and 𝑑 = 10 𝑚𝑚.) The convergence can usually be reached at iteration numbers less than 100. The entire simulation time on a regular modern PC is around 1-3 seconds when this model is implemented in MATLAB. The 𝑛 = 1 case corresponds to the Laplace equation's solution. Initially, the space charge motive is underestimated and therefore leads to high electron density. The exaggerated density causes higher motive due to space charge (𝑛 = 6) which in turn reduces the electron density in the subsequent iterations until convergence. 49 Figure 3-2. The evolution of the motive vs. distance as the loop in Figure 3-1 is repeated. 𝑛 represents the iteration number. The 𝑛 = 1 case corresponds to the Laplace equation solution. This initial attempt underestimates the space charge and therefore leads to high electron density. The exaggerated density results in higher space charge effect (𝑛 = 6), which in turn reduces the electron density in subsequent iterations. The parameters of the TEC under study are 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 2,000 𝐾 and 𝑑 = 10 𝑚𝑚. 3.2.1. Particle tracing approach to bypass the analytical solution to the Vlasov equation Here, we show that particle tracing can be used to bypass the analytical solution of the Vlasov equation. A finite difference time-domain (FDTD) particle-in-cell (PIC) approach to obtaining the output characteristics of thermionic converters has been adopted by Lo et al. using a 1-dimensional object oriented particle tracing software package called OOPD1[65]. This method can be computationally expensive, since the electrons are followed in real time and their influence on each other is accounted for by solving the Poisson equation. Since we are interested in the steady state response of the system, we propose a particle tracing approach that can mimic the solving of the Vlasov equation and integration of the distribution function in an efficient manner. Instead of solving the Poisson equation at each time step as the electrons propagate, we solve the Poisson equation at the end of each iteration where the equilibrium electron density 50 for that particular motive profile has been reached. Therefore, in this approach, we are integrating the equations of motion as a way to find solutions to the Vlasov equation and use time averages of the desired variables (such as electron density) to represent the values of these variables that are obtained from the actual distribution functions. The premise of this approximation is based on the ergodicity of the phase space of the electrons [50]. Another approach could be to solve the Vlasov equation as a highly non-linear partial differential equation. Though this strategy is feasible in 1-D, leading to two variables, it can quickly become intractable in 3-D cases where the distribution function lives in the vast 6-D space. The difficulty of the particle tracing approach on the other hand scales linearly with the dimensionality of the phase-space. Therefore, adopting a particle tracing approach is practical since it is possible to sample the phase space more easily by following a number of test particles in each step. A flux of 107 electrons is released from the emitter at regular intervals, over a total duration of 5 𝜇𝑠 (roughly 1,000 times larger than the time-of-flight of the electron with an average velocity corresponding to 2,000 𝐾). The initial velocities of these electrons follow an HM distribution. The effective area of emission is adjusted so that the flux of the forward-moving electrons at the emitter matches the incident thermionic flux, 𝐽𝑠𝑎𝑡/𝑒 . Using a higher number of incident electrons leads to a higher number of electrons in the inter-electrode region, although the overall electron density remains constant. These electrons are followed on their trajectory to the collector through the potential landscape produced in a previous step as the solution of the Poisson equation. The acting force at each particle position is calculated at a grid point and interpolated to the current position of the particle. The particle is propagated using the velocity Verlet algorithm [89] [90] [50]: 𝑥𝑖 = 𝑥𝑖−1 + 𝑣𝑖−1 Δ𝑡 +12𝐹(𝑥𝑖−1)𝑚Δ𝑡2 and 𝑣𝑖 = 𝑣𝑖−1 +12(𝐹(𝑥𝑖−1)+𝐹(𝑥𝑖))𝑚Δ𝑡, where 𝑥𝑖, 𝑣𝑖 and 𝐹(𝑥𝑖) are the position, velocity and the acting force on the electron, respectively at step 𝑖 and Δ𝑡 is the time interval. Each electron “feels” the presence of the other electrons only though the Poisson equation in a mean field fashion. The electrons that are pushed back to the emitter or reach the collector are removed from the system and their numbers are recorded (along with their arrival time) and used subsequently to calculate the current density. If the final position of an electron lies between two grid points, its charge is distributed between its nearest neighbors using a linear weighting function. The electron density obtained by coupling the particle tracing approach and the Poisson equation were compared with that of Vlasov-Poisson system. The results are depicted in Figure 3-3(a), when the solution has been converged after 20 rounds of iteration. The fluctuations in the electron density are due to the element of randomness in the initial velocity of the electrons introduced by the random Gaussian distribution that was used to generate the initial velocities. However, these fluctuations are washed out in the double integration process of solving the Poisson equation, leading to the same motive and hence the same current (Figure 3-3(b)). The number of the emitted electrons should be chosen judiciously to ensure that the number of electrons in the inter-51 electrode region is higher than the number of grid points; otherwise, the resulting electron density becomes disjointed. Figure 3-3. Comparison between the particle tracing approach to calculate the electron density and the analytical Vlasov solution. The average electron density (a) and motive (b) for the final round of the iteration, when the convergence is reached, are plotted here. In the particle tracing approach, 107 particles are tracked in each iteration. This plot corresponds to the converged solution after 20 rounds of iteration. The electron density and motive plots from particle tracing are calculated based on their average steady-state values. (c) Comparison between the improved-and-extended-Langmuir (IEL) and the numerical 52 solutions (i.e. using Vlasov integrals or particle tracing) for the current-voltage characteristics of the same TEC as in part a) with 𝜙𝐸 =4.5 𝑒𝑉, 𝑇𝐸 = 2000 𝐾 and 𝑑 = 1 𝑚𝑚 (The overall current is calculated for an emission surface area of 3.14 × 10−8 𝑚2). The IEL results were obtained based on the strategy outlined in [11]. The particle tracing current is the average steady state and matches the solutions obtained by applying an asymptotic expansion to the electron density. It is noted that this approach is more time-consuming (about 2-3 hours on a regular PC) than using an analytical solution to the Vlasov equation (several seconds). However, the particle tracing approach is not limited to cases where an analytical solution to the Vlasov equation exists. This model was applied to a wide range of the parameters of TECs, including different emitter temperatures, inter-electrode distances and workfunctions. We used the analytical solutions that exist in this special case based on the IEL approach that we proposed in [2] to measure the numerical accuracy of this new self-consistent method. As depicted in Figure 3-3(b), the results from this method are in complete agreement with those of the IEL approach. The parameters of the TEC under study are 𝜙𝐸 = 𝜙𝐶 =4.5 𝑒𝑉, 𝑇𝐸 = 2,000 𝐾 and 𝑑 = 1 𝑚𝑚. Furthermore, the solutions obtained by the asymptotic expansion of the electron density obtained from the Vlasov equation and the particle tracing approach completely agree with each other. Self-consistent solution in the presence of back-emission This model can be naturally extended to include the effect of electron emission from the collector as well. It is assumed that some electrons are emitted from the collector and are fully absorbed as they reach the emitter. The interaction between the electrons ejected from the emitter and collector in the inter-electrode region occurs through the Poisson equation; the motive is cast in the form 𝜓(𝑥) = 𝜓𝑙(𝑥) +𝜓𝑑𝑝,𝐸(𝑥) + 𝜓𝑑𝑝,𝐶(𝑥), where the last two terms represent the solutions to the Poisson equation with homogeneous Dirichlet boundary conditions for electron densities originating from the emitter and collector, respectively. The ranges of velocities of the electrons originating from the collector can be worked out by using arguments similar to those in the previous section. Recalling that positive velocity is defined along the positive direction of the x axis, the maximum velocity is positive for electrons at 𝑥 ≥ 𝑥𝑀 and equal to (2 𝜓𝑀−𝜓(𝑥)𝑚)12, which is the same as 𝑣𝑥,𝑚𝑖𝑛 for electrons originating from the emitter; electrons originating from the collector can span the range from −∞ to 𝑣𝑥,𝑚𝑖𝑛. On the other hand, at 𝑥 < 𝑥𝑀, the maximum 53 velocity of the electrons originating from the collector is negative and equal to −𝑣𝑥,𝑚𝑖𝑛 = −(2 𝜓𝑀−𝜓(𝑥)𝑚)12. Therefore, the overall electron velocity distribution can be written as 𝑓(𝑥, 𝑣) = 2 𝑛𝐸(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) ×Θ(𝑣𝑥 ∓ 𝑣𝑥,𝑚𝑖𝑛)+2 𝑛𝐶(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐶)32exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐶−𝑚𝑣22 𝑘𝐵𝑇𝐶) ×Θ(−(𝑣𝑥 ∓ 𝑣𝑥,𝑚𝑖𝑛)), (3-5) for 𝑥 ≥ 𝑥𝑀 (top sign) and 𝑥 < 𝑥𝑀 (bottom). 𝑛𝐸(𝑥𝑀) and 𝑛𝐶(𝑥𝑀) represent the contributions of electrons arising from the emitter and the collector, respectively, to the electron density at the point of maximum motive. Equation (3-5) is subsequently integrated to obtain the electron densities and replaced in the Poisson equation, 𝑑2𝜓(𝑥)𝑑𝑥2= −𝑒2𝜖0𝑛(𝑥𝑀) {exp(𝛾𝐸) (1 − erf (𝛾𝐸12)) + exp(𝛾𝐶) (1 + erf (𝛾𝐶12)) ,exp(𝛾𝐸) (1 + erf (𝛾𝐸12)) + exp(𝛾𝐶) (1 − erf (𝛾𝐶12)) ,𝑥 > 𝑥𝑀𝑥 ≤ 𝑥𝑀, (3-6) where 𝑛(𝑥𝑀) = 𝑛𝐸(𝑥𝑀) + 𝑛𝐶(𝑥𝑀) and 𝛾𝐶 ≡ (𝜓𝑀 − 𝜓(𝑥))/𝑘𝐵𝑇𝐶. This equation is solved by the same strategy as in the last section. Additionally, by arguments similar to those preceding equation (3-4), the overall current density at each iteration number is calculated as 𝐽 = 𝐽𝑠𝑎𝑡,𝐸 exp (−𝜓𝑀𝑘𝐵 𝑇𝐸) −𝐽𝑠𝑎𝑡,𝐶 exp (−𝜓𝑀+𝑒 𝑉𝑘𝐵 𝑇𝐸), where 𝑉 is the potential difference between the points just outside the collector and emitter. The first term represents the contribution of the emitter to the total current density, 𝐽𝐸, whereas the second term is the contribution of the collector, 𝐽𝐶. 54 55 56 Figure 3-4. (a) The contributions of the emitter and collector current densities, 𝐽𝐸 and 𝐽𝐶, to the total current density, 𝐽𝑡𝑜𝑡 of a TEC with 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 2,000 𝐾, 𝑇𝐶 = 1,900 𝐾 and 𝑑 = 1 𝑐𝑚. (b) The evolution of the total current density as a function of the iteration number, 𝑛. (c) The converged motive as a function of position for an applied voltage of 1 𝑉. (d) The contributions of the electrons originating from the emitter and the collector to the overall electron density as a function of position. The value of 𝑛(𝑥𝑀) is calculated by individually equating 𝐽𝐸 and 𝐽𝐶 to the integral of the x component of the velocity and the velocity distribution functions, equation (3-5), resulting in 𝐽𝐸 =𝑒 𝑛𝐸(𝑥𝑀) (2 𝑘𝐵 𝑇𝐸𝜋𝑚)1/2and 𝐽𝐶 = 𝑒 𝑛𝐶(𝑥𝑀) (2 𝑘𝐵 𝑇𝐶𝜋𝑚)1/2. In the last two sections, it was assumed that electrons maintain their temperatures as they propagate to the opposing electrodes. Since 𝐽𝐶 and 𝐽𝐸 depend on 𝑇𝐸 and 𝑇𝐶, and the equation of continuity of charge dictates that 𝛁. 𝑱 be constant in the steady-state, it follows that the electron temperature remains constant. The validity of this argument is nonetheless dependent on how well the velocity distribution of the electrons can be captured by an HM distribution at the point of the maximum motive. To calculate a limit on the validity of the HM distribution, the differential scattering cross-section formula [91] for the Coulomb interaction potential between two electrons was considered. The total scattering cross-section, 𝜎𝑡𝑜𝑡, can be calculated by integrating the differential scattering cross-section, resulting in infinity. However, the HM approximation becomes invalid when the scattering angle is more than 90 °. Using this angle as the lower bound leads to the following non-HM total scattering cross-section: 𝜎𝑡𝑜𝑡,𝑛𝐻𝑀 = 𝜋 (𝑒22 𝜋 𝜖0 𝑚 𝑣𝑟𝑒𝑙2 )2, where 𝑣𝑟𝑒𝑙 is the relative velocity between the electrons. The distance that electrons travel before this effect is considerable can be approximated as (𝑛(𝑥) 𝜎𝑡𝑜𝑡,𝑛𝐻𝑀)−1 . In the case of the device examples analyzed in this chapter, this distance is ~ 700 𝜇𝑚. Even for a high-performance TEC with 𝜙𝐸 = 2.0 𝑒𝑉 and 𝑇𝐸 =1600 𝐾, this distance is ~ 200 𝜇𝑚. Given that a modern TEC can easily have an inter-electrode distance smaller than this, our model has a very broad range of applicability. The contributions of the emitter and collector current densities, 𝐽𝐸 and 𝐽𝐶, to the total current density, 𝐽𝑡𝑜𝑡, of a TEC with 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 2000 𝐾, 𝑇𝐶 = 1900 𝐾 and 𝑑 = 1 𝑐𝑚 are presented in Figure 3-4 (a). The evolution of the total current density as a function of the iteration number, 𝑛, is plotted in Figure 3-4 (b) In the case of an applied voltage of 1 𝑉, the converged motive as a function of position is displayed in Figure 3-4 (c). Finally, the contributions of the electrons arising from the emitter and the collector to the overall electron density are depicted in Figure 3-4 (d). Expectedly, the electron density arising from the emitter has its maximum around the emitter, whereas the electrons originating from collector are mostly concentrated in the proximity of the collector. Therefore, the overall electron density has two maxima in the presence of back-emission. 57 Lastly, it is noted that the methods developed in this Chapter can be used to calculate the output characteristics of TECs for a wide range of workfunctions, temperatures, applied voltages and inter-electrode distances (Figure 3-5). This model applies as long as the quantum effects (e.g. tunneling in nanoscale gaps) and relativistic effects (at extremely high biases and not applicable in the power generation mode) are negligible. The cases studied here were deliberately chosen to have high space charge (and naturally small output power density) so that electron transport in these conditions could be investigated. The maximum output power density as a function of inter-electrode gap size is shown in Figure 3-6. It is assumed that 𝜙𝐸 = 2.5 𝑒𝑉, 𝜙𝐶 = 2.0 𝑒𝑉, 𝑇𝐸 = 1,800 𝐾 and 𝑇𝐶 = 700 𝐾. Based on the proposed method, for each value of gap size, the output current density for various applied voltages was calculated and subsequently used to calculate the power density. The maximum value of power density at each gap size (which occurs at different voltages for different gaps) was used to produce the data in Figure 3-6. Expectedly, at smaller gaps where the space charge effects are less palpable, the change in the maximum power density as a function of distance is slower than at the higher distances, where the space charge effects are more prominent. These trends are in accordance with the power densities calculated in other works [14]. 58 59 Figure 3-5. The characteristics of a thermionic converter in the presence back-emission. The workfunctions of the emitter and the collector are both 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, the emitter’s temperature is, 𝑇𝐸 = 1500 𝐾, the collector’s temperature is 𝑇𝐶 = 1400 𝐾 and the inter-electrode distance is 𝑑 = 100 𝜇𝑚. Evolutions of the (a) current density from the emitter, 𝐽𝐸, (b) motive, 𝜓 (c) electron density, 𝑛𝐸,𝑒(𝑥) from the emitter and (d) electron density from the collector, 𝑛𝐶,𝑒(𝑥), are plotted. In part (a), 𝛿 is the mixing ratio between motive in subsequent iterations of solving the Vlasov-Poisson system 60 Figure 3-6. The maximum output power density of a TEC as a function of the inter-electrode distance. The device parameters are 𝜙𝐸 = 2.5 𝑒𝑉, 𝜙𝐶 = 2.0 𝑒𝑉, 𝑇𝐸 = 1,800 𝐾 and 𝑇𝐶 = 700 𝐾. The proposed method was employed to calculate the output current density at different voltages. These values were subsequently used to calculate the power density. The maximum power density was found by sweeping the voltage and calculating the current at each value of the gap size, due to the space charge effect. The values calculated in this figure are in accordance with the values reported in the literature [14], [57]. Convergence of the solutions to the Vlasov-Poisson system As with any two coupled differential equations with necessarily shared parameters, one expects that three scenarios could arise: the solutions to the two equations diverge, oscillate or converge in each subsequent step. Divergence, in an absolutely mathematical sense, could occur if the two equations are incompatible with each other, i.e., the increments at each step amplify each other because a common solution to the two equations does not exist. Another type of diverging solution can arise even for compatible equations when the damping ratio is not chosen properly and the intermediary solutions exceed 61 the maximum finite floating-point number permissible in the programming language. The other possibility is for the solutions to oscillate, which could occur in cases where a common solution to the two equations does not exist, although the changes do not amplify each other (i.e., the electron density is not an increasing function of motive, whereas the motive is an increasing function of electron density), or if a small value of a damping ratio is chosen. Lastly, the equations can converge if they share a common solution and if a proper damping ratio is selected. The diverging solutions should be ruled out in the Vlasov-Poisson system, since the changes in one do not amplify the solutions of the other; e.g., incrementing the motive, would significantly reduce the number of electrons with enough kinetic energy to overcome the maximum barrier, leading to reduced potential in the next step. Therefore, it remains to investigate if the Vlasov-Poisson system belongs to the oscillating type or the converging type. It can be seen that the steady state response of a TEC in the absence of ions or alternating potentials (all the cases studied in this Chapter) belongs to the converging category by performing the following analysis: 1. Choose a value of 0 < α < 0.1 and run the simulation. 2. If the solution converged after a certain iteration number, n, set α to a higher value, i.e., 0.2 and run the simulation again to a high number of iterations (e.g., n = 1000). If the solution starts to oscillate, it shows that the system does not have a common solution to both equations. However, if the solution remains constant in the succeeding iterations, the system has settled on the common solution between the two equations. The latter occurs in all cases that we have studied. The mixing ratio, α, cannot be set to its maximum, i.e., 1. This is because the solution to the Poisson equation is calculated numerically and its accuracy is determined by the fineness of the mesh. The slightest numerical mismatch in the motive (an error on the order of 10−10) can grow exponentially due to the highly nonlinear nature of the Vlasov equation. For our cases, the number of mesh elements in the inter-electrode region is 1000. A mixing ratio of 0.2 was found to be enough to dampen the numerical errors that arise from the finite size of the mesh elements plus all the numerical errors in the other calculations. 3. If the solution oscillates and never reaches convergence, use a smaller value of α, but continue the simulation to a higher number of iterations, since the changes at each step have been damped more severely. 4. Go back to 2 (if the solution converged) and 3 (if the solution was oscillatory). If the solution is always oscillatory (step 3), it can be imagined that a common solution to the two equations does not exist. This has never been observed for the cases that we have studied in this chapter. A complete proof of the oscillatory nature of the solution would require an infinite number of iterations, with 62 ever smaller values of α. This is obviously not feasible using a computer since the smallest floating-point number of the programming language would soon be reached. However, this has never been observed in our cases (depending on the value of α, the test ends immediately in step 2 with non-oscillatory solutions, or for higher values of α in step 3 for one iteration and subsequently in step 2 with a non-oscillatory solution.) Therefore, this could be considered an a posteriori proof for the convergence of the solutions of the Vlasov-Poisson system in the steady state response of the TECs in the absence of ions or alternating currents. An a priori argument for the existence of converging solutions to the Vlasov-Poisson system in the space charge limited regime can be deduced from the work of Langmuir and Hatsopoulos [56] [57], where they present a closed-form solution to the system, signifying a unique current associated with each voltage in the space charge limited mode, in the absence of the back-emission. In the saturation regime, the current is constant (for as long as the Schottky effect is negligible). In the retarding mode, there is also a unique correspondence between the current and the applied voltage. The 1-D model can in some cases, capture cases other than two infinitely large parallel electrodes, such as the case where spherical or cylindrical symmetries exist. However, in more general cases, the model should be extended (as indicated in the future works section), to include extra dimensions which allows calculation of the angular distribution function of the emitted electrons. Summary We presented a self-consistent approach to derive the characteristics of thermionic energy converters. In the case that the analytical solution of the Vlasov equation is not possible, our model uses a particle tracing approach to calculate the electron densities. Moreover, by introducing asymptotic expansion to deal with the ill-behaving functions, the range of applicability of our approach is significantly increased in the presence of an analytical solution. The model also circumvents the splicing issues and captures the entire characteristics of thermionic converters. The results of this strategy were shown to be in agreement with previous solutions in the absence of back-emission. Subsequently, this method was employed to calculate the output characteristics of thermionic converters in the presence of back-emission. 63 Chapter 4- Analysis of thermionic energy converters comprising a grid4 Quantum mechanics tells us that the bound states of a potential well are quantized—a phenomenon that is easily understandable based on wave properties and resonance. Here, we demonstrate a classical mechanism for the formation of a momentum gap in the phase-space of electrons traveling as particles in a potential well in vacuum. This effect is caused by the reflection of electrons from at least two potential maxima, which might, for instance, exist due to space charge distribution in a triode configuration. This gap plays a critical role in space charge-mitigated electron transport in vacuum, such as in a thermionic energy converter with a positively biased grid (Figure 1-6), where we show that the current density can be increased by 3 orders of magnitude. Theory and model In this work, it is assumed that the grid pitch is finer than the Debye length of the electron gas (typically 10 − 500 𝜇𝑚). Under this condition, we assume that the plane of the grid (including the holes) can be considered to be approximately an equipotential surface and thus the motive at the gate position can be set equal to −𝑒 𝑉𝐺, where 𝑉𝐺 is the grid bias with respect to the emitter. At first, we will develop the solutions to the Vlasov equation for different operation regimes and potential schemes in the presence of the grid. These velocity distributions are then employed to calculate the electron density and the motive profile by solving the Poisson equation. Vlasov equation in the presence of a grid At first we will provide a brief introduction about the behaviour of a TEC in the case that a grid is present and under such conditions that a momentum gap arises in the phase-space of electrons. This introduction is followed by a detailed analysis of all the different motive scenarios that could arise in the presence of a grid. 4 This chapter is reproduced with minor modifications from [4]. 64 Figure 4-1. (a) The potential (𝜓) or motive profile of a potential well as a function of distance, 𝑥 and the resulting velocity distribution functions, 𝑓(𝑥, 𝑣𝑥) in different regions. (b) Schematic of a vacuum device comprising of an emitter, a collector and an ancillary electrode (grid) that allows the transport of electrons through the classical potential well. The potential well depicted in (a) can be generated due to the external applied fields as well as the interaction between the electrons in the inter-electrode region. 65 Imagine that a stream of electrons with a wide range of positive initial velocities is fired (from the left hand side) at the potential (motive) profile depicted in Figure 4-1. (a) and is fully collected once it reaches the right hand side. For the moment, the interactions among the electrons are neglected; they will be subsequently included in the more detailed discussion below. On their path, electrons face two potential barriers in the inter-electrode region, and a potential well between those. If the height of the second barrier is higher than that of the first one, the reflection of the electrons from these two barriers will lead to an interesting phenomenon, whereby certain momenta will not be permissible in the potential well between points 2 and 5. The emergence of this behaviour can be understood by considering that only electrons with initial kinetic energies above the first barrier (point 2) can enter the well. At any subsequent point in the well, these forward-moving electrons will have velocities greater than a certain value v2(x), which depends on the local potential. These constitute the positive, semi-infinite shaded band on the corresponding velocity distribution graph on Figure 4-1. (a). Of these electrons, those with a kinetic energy below the second barrier will be reflected and constitute the negative, finite shaded band on the velocity distribution graph. This leads to a position-dependent momentum gap centered around zero, as can be seen on the distribution graph. In the regions between points 1 and 2, and between points 5 and 6, which are outside the potential well, this gap disappears as the electrons, both forward-moving and those reflected from the barriers, can have velocities all the way down to zero. Similarly, in the region between points 6 and 7, electrons will have only forward momentum and there will be no gap. Note that, if the first potential barrier (point 2) is taller than the second one (point 6), this phenomenon will not occur, since in this case only electrons with initial kinetic energies higher than the tallest barrier will enter the potential well and they will not be reflected from the second barrier. It is interesting to observe that the momentum gap discussed is a result of having forward-moving and reflected electrons in the system. This is reminiscent of the quantization of momentum in a quantum well, which could be regarded as the result of interference of forward-travelling and backward-travelling waves. Nevertheless, the two are fundamentally different effects. We now turn to a detailed analysis of the problem and discuss the profound role of this momentum gap on electron transport. The potential well portrayed in Figure 4-1. (a) can be generated by the device configuration depicted in Figure 4-1. (b). Electrons with various kinetic energies are released from a source electrode (or emitter), pass through a potential well and are collected at another electrode (the collector). The potential well can be formed by the application of one or several auxiliary electrodes, namely, grids (Figure 4-1. (b)), or the motive can be lowered by other means, such as the presence of positive ions in the system (Chapter 5), or the use of negative electron affinity material [92]. An immediate way for realizing 66 different initial velocities at the emitter surface is the ejection of electrons from a hot electrode based on thermionic emission. Such electrons assume a hemi-Maxwellian (HM) distribution in steady state at the emitter if the potential of the point outside the collector lies below that just outside the emitter. Here, we show how these effects can be quantified by solving the Vlasov equation. The electrodes are assumed to be infinitely large in lateral dimensions, so that a 1-D Vlasov equation can be used. (We expect that the momentum gap will leave its mark in the 3-D case as well, although the analysis would be more complex in that case.) The virtue of using the Vlasov equation is that the electron-electron interactions are included (albeit in a mean-field manner), leading to correct results in the steady state when coupled to the Poisson equation. To keep the notation uncluttered, some new variables are introduced in this section to provide the solutions to the Vlasov equation (1-3).The solution to this equation can assume the form of any arbitrary function of an integral of motion [51], such as the energy of an electron, 12𝑚𝑣2 + 𝜓(𝑥). The exact solution can be obtained by applying the proper boundary condition (Figure 4-1(a)), i.e., electrons having an HM distribution at the point of maximum motive, (𝑥𝑀 , 𝜓𝑀), which can be written as 𝑓(𝑥𝑀 , 𝑣) =2 𝑛(𝑥𝑀)𝑣𝑡ℎ−3 exp (−𝛽 (12𝑚𝑣2 + 𝜓(𝑥𝑀) − 𝜓𝑀)) × Θ(𝑣𝑥), where 𝜓(𝑥𝑀) ≡ 𝜓𝑀, 𝑣𝑡ℎ is the average thermal velocity of the HM distribution, (2 𝜋𝑚𝛽)12, Θ is the unit step function, and 𝛽 =1𝑘𝐵𝑇𝐸 (cf. equation 1-5). Depending on the actual voltage and current values in the device, a total of 8 possible motive configurations leading to different distribution functions could arise in the presence of an emitter, a collector and an ancillary electrode. A case that satisfies the two aforementioned conditions required for the emergence of velocity gaps is depicted in Figure 4-1.(a). Electrons on the right hand side of the maximum motive (points 6-7) can only move unidirectionally since they had enough kinetic energy to surmount the motive barrier, with a minimum velocity 𝑣𝑚𝑖𝑛 = ( 2𝛾𝛽𝑚)12, where 𝛾 ≡ 𝛽(𝜓𝑀 − 𝜓(𝑥)). Therefore, for such points, 𝑓(𝑥, 𝑣) =2 𝑛(𝑥𝑀)𝑣𝑡ℎ−3 exp (−𝛽 (12𝑚𝑣2 + 𝜓(𝑥) − 𝜓𝑀)) × Θ(𝑣𝑥 − 𝑣𝑚𝑖𝑛). After integration, this results in the electron density 𝑛(𝑥) = 𝑛(𝑥𝑀) exp(𝛾) (1 − erf 𝛾1/2). In regions where electrons face a barrier (or barriers) only in the forward direction, that is between points 1 and 2, and between points 5 and 6, they move bidirectionally since those without sufficient energy to overcome all the barriers are scattered and fully reverse their direction, resulting in negative velocities. In such cases we have 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀)𝑣𝑡ℎ−3 exp (−𝛽 (12𝑚𝑣2 + 𝜓(𝑥) − 𝜓𝑀)) × Θ(𝑣𝑥 + 𝑣𝑚𝑖𝑛) and 𝑛(𝑥) = 𝑛(𝑥𝑀) exp(𝛾) (1 + erf 𝛾1/2). A special behaviour arises at points between the two peaks 𝜓𝑀,𝑙 and 67 𝜓𝑀, with motives smaller than 𝜓𝑀,𝑙 (points 2-5). In such cases, the first, smaller peak blocks the electrons with initial kinetic energies smaller than 𝜓𝑀,𝑙. This results in 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀)𝑣𝑡ℎ−3 exp (−𝛽 (12𝑚𝑣2 +𝜓(𝑥) − 𝜓𝑀)) × { Θ(𝑣𝑥 − 𝑣𝑚𝑖𝑛) − Θ(𝑣𝑥 + 𝑣𝑚𝑖𝑛) + Θ(𝑣𝑥 + 𝑣2)}, where 𝑣2 = ( 2𝛾𝑙𝛽𝑚)12and 𝛾𝑙 ≡ 𝛽 (𝜓𝑀,𝑙 −𝜓(𝑥)), according to Figure 4-1. (a). The velocity gap centered around zero is due to the absence of electrons with initial kinetic energies between 0 and 𝜓𝑀,𝑙, whereas the absence of velocities less negative than -𝑣𝑚𝑖𝑛 is because any electron with initial kinetic energy higher than the tallest peak (𝜓𝑀,𝑙) will not be reflected. This leads to the distinct, negative velocity band seen on Figure 4-1. (a) in the region between points 2-5. Note that this band structure depends on position and the gap is widest at the location of the grid electrode. Integration of the distribution function results in 𝑛(𝑥) = 𝑛(𝑥𝑀) exp(𝛾) (1 − 2 erf 𝛾𝑙1/2+ erf 𝛾1/2). Electron densities at high values of 𝛾 were calculated by using an asymptotic expansion of the error function [88] or numerically integrating the velocity distribution function using the method of global adaptive quadrature [93]. Other possibilities of the motive profiles are considered now. Typical possibilities of the potential profile that could arise as a grid is inserted in the inter-electrode distance of a TEC are plotted in Figure 4-2.(b-h). In Figure 4-2.(b) the motive of the point just outside the collector lies below the motive of the point just outside the emitter. In other words, the internal voltage difference, 𝑉𝑖𝑛𝑡, between the collector and the emitter is positive. In the case that the workfunction of the emitter and the collector are equal, the internal voltage difference is equivalent to the applied bias, 𝑒 𝑉. In the more general case, according to Figure 4-2, 𝑉𝑖𝑛𝑡 = 𝑉 + (𝜙𝐸 − 𝜙𝐶)/𝑒. Moreover, the space charge effects due to the presence of electrons in the inter-electrode region are not strong enough to produce a maximum motive in the inter-electrode region, but merely curve the potential profile upwards (due to the negative charge of the electrons). In this case, all electrons with enough energy to circumvent the emitter’s workfunction can make it to the other side, i.e., undergo a unidirectional transport in the inter-electrode region. Therefore, according the Figure 4-2.(b), the minimum velocity of the electrons at an arbitrary point, x, is 𝑣1 = (2 −𝜓(𝑥)𝑚)12 and the velocity distribution is equal to 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥1) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 − 𝑣1). 𝑛(𝑥1) is the electron density at point 1 (𝑥 = 0). Expectedly, at 𝑥 = 0, 𝑣1 = 0 and the probability distribution assumes an HM distribution. Figure 4-2.(c) depicts a similar potential profile, but with drastically different electron statistics. This behaviour is significantly different from the statistics in parts a and b of Figure 4-2. This distribution can 68 be understood by considering the initial kinetic energies of the electrons at the point just outside the emitter (point 1). At an arbitrary position, x (such a point 2 or 3), between points 1 and 4, the electrons can have both positive and negative velocities. The starting point of the trajectory of the electrons is at point 1, and the minimum kinetic energy is equal to zero. These electrons, with zero kinetic energy will have a velocity equal to 𝑣1 = (2 −𝜓(𝑥)𝑚)12 at point x. Subsequently, these electrons can travel as far as point 4, where their kinetic energy becomes zero. These electrons are then reflected by the motive at points following 4. The velocity of these electrons, on their return journey towards the emitter, is equal to −𝑣1. Electrons at point 1 with velocities slightly higher than zero will have a velocity higher than 𝑣1 at point x, since their initial kinetic energy was higher. Additionally, these same electrons, after reflection, will have a more negative velocity on their way back. This means that there is a velocity gap at point x for velocities, 𝑣𝑥, that satisfy −𝑣1 ≤ 𝑣𝑥 ≤ 𝑣1. This velocity gap is not occupied by electrons and needs to be carved out of the Maxwellian distribution. The other part of the velocity distribution in Figure 4-2.(c) that needs to be excluded can be understood as follows: at point 1, electrons that have kinetic energies higher than the point of maximum motive, i.e., point 6 with a motive equal to −𝑒 𝑉𝑖𝑛𝑡 will not be reflected and can travel unidirectionally to the collector. This corresponds to a velocity, 𝑣2 = (2 −eV𝑖𝑛𝑡−𝜓(𝑥)𝑚)12 at point x. Electrons with energies less than this threshold travel bidirectionally. Therefore, an additional velocity gap from −∞ to −𝑣2 needs to be excluded from the Maxwellian velocity distribution. Therefore, the overall velocity distribution function can be written as 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥6) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (−eV𝑖𝑛𝑡−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) ×{ Θ(𝑣𝑥 − 𝑣1) − Θ(𝑣𝑥 + 𝑣1) + Θ(𝑣𝑥 + 𝑣2)}. For points 4-6, the motive is positive and therefore all the incoming electron energies from 0 to −𝑒 𝑉𝑖𝑛𝑡 can fill the velocity distribution function. The presence of the middle gap in the velocity distribution function for points 1-4 is due to the fact that no incoming electron energies correspond to negative values of the motive. Obviously, this phenomenon does not arise in positive motives. Therefore, for points 4-6, the velocity distribution function is given as 𝑓(𝑥, 𝑣) =2 𝑛(𝑥6) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (−eV𝑖𝑛𝑡−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 + 𝑣2). As a check to validity of this formulation, it is noted that these two equations (governing the regions between points 1-4 and 4-6, respectively) yield the same velocity distribution at their common point, i.e., point 4. Additionally, the value of 𝑣2 at point 6 is zero, leading to the familiar HM distribution at this point. Also, it is noted that the points 2 and 3, by having the exact same motive, are equivalent with each other and have the same velocity distribution. This does not apply in the case of points 2 and 4 in Figure 4-2.(a), where although the points have the same motive, their electrons follow different statistics. 69 70 71 72 Figure 4-2. (a) The motive diagram of a regular TEC in the space charge limited mode. The possible velocity distribution functions are plotted as well. The electrons can maintain different velocity distributions depending on which side of the 𝜓𝑀 they are positioned at. The possible motive potentials and their corresponding velocity distribution functions are depicted in (b)-(h). (b) The point of maximum motive coincides with the point right outside the emitter. In this case, the electrons can only move unidirectionally. (c) The point of maximum motive occurs at the collector, and the space charge is not strongly felt. The possible velocity distribution functions are plotted underneath the motive profile. d) Space charge is more intense in the emitter-grid region. Electrons with enough energy to circumvent the maximum motive can move only unidirectionally. (e) The point of maximum motive resides on the collector, though the space charge is still more severely felt in the emitter-grid region. (f) The space charge is more intense on the grid-collector side and the maximum motive lies in this region. (g) The space charge is tangible in both emitter-grid and grid-collector regions, with the motive height at its maximum peak on the grid-collector side. (h) Same as g, save that the maximum peak occurs on the emitter-grid side and the peak on the grid-collector side serves as a local extremum point. 73 In Figure 4-2. (d), the amount of space charge is high enough to lead to a peak in the distance between the emitter and the grid. The space charge electron density only curves up the motive potential in the region between the grid and the anode without producing a peak. Moreover, the value of the motive at point 3 is higher than the motive at point 10 due to the internal motive difference, −𝑒 𝑉𝑖𝑛𝑡. This configuration follows a similar statistics as a regular TEC, since the peak occurs before the grid. Therefore, electrons that have enough kinetic energy to overcome the maximum motive are not reflected and traverse unidirectionally towards the collector. For points 1-3, the electrons move bidirectionally and their distribution can be obtained from the positive sign in the argument of the unit step function in the distribution function, and 𝑣1 = (2 𝜓𝑀−𝜓(𝑥)𝑚)12. For points 3-10, the distribution function is given by the negative sign of the argument of the unit step function with the same value of 𝑣1. The motive profile in Figure 4-2. (e) is similar to part d, save that the maximum motive lies on the point just outside the collector (point 10). Thus, the extremum point 3 is a local maximum with coordinates 𝑥𝑀,𝑙 and 𝜓𝑀,𝑙. This again leads to a substantially different electron statistics than that of Figure 4-2. (d) In the region between points 1-3, the electrons traverse bidirectionally with their velocity distribution given by 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥10) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (−eV𝑖𝑛𝑡−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 + 𝑣1), and 𝑣1 is (2 –eV𝑖𝑛𝑡−𝜓(𝑥)𝑚)12. At points 3-9, though the electrons have passed the local barrier at 𝑥𝑀,𝑙, they still traverse bidirectionally, since their energy may only be sufficient to overcome the barrier 𝜓𝑀,𝑙, but not high enough to overcome the motive at point 10, −𝑒 𝑉𝑖𝑛𝑡. The minimum kinetic energy that the electrons need to possess at the point just outside the emitter, point 1, in order to reach to point 3, is equal to 𝜓𝑀,𝑙. Once these electrons are at an arbitrary point x between points 3-9 their positive velocity is 𝑣2 = (2 𝜓𝑀,𝑙−𝜓(𝑥)𝑚)12, when they initially pass across point x. However, these electrons do not have enough energy to reach point 10. These electrons are thus reflected towards the emitter and reach point x with a negative velocity equal to −𝑣2. Electrons whose absolute value of their velocities are less than 𝑣2 cannot be present at point x. Therefore, a slice of – 𝑣2 ≤ 𝑣𝑥 ≤ 𝑣2 is removed from the Maxwellian distribution. Another part that needs to be excluded from the Maxwellian distribution is the contribution of the electrons with velocities more negative than −𝑣1, since these electrons already possess enough kinetic energy to overcome the point of the maximum motive at 10. Therefore, the overall velocity distribution can be written as 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥10) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (−eV𝑖𝑛𝑡−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × { Θ(𝑣𝑥 − 𝑣2) − Θ(𝑣𝑥 + 𝑣2) + Θ(𝑣𝑥 +𝑣1)}. 74 If electrons have made it to the region between points 9 and 10, their initial kinetic energy is more than 𝜓𝑀,𝑙 and the electrons are allowed to have an initial kinetic energy that is equal to 𝜓(𝑥) and thus the middle gap (– 𝑣2 ≤ 𝑣𝑥 ≤ 𝑣2 ) in the Maxwellian distribution disappears. Therefore, in the region 9-10, the electrons follow the same statistics as that of points 1-3. Importantly, the velocity distribution function at the boundary points 3 and 9 yield the same results from the two different distributions. Moreover, at point 10, the electrons form an HM distribution. The maximum motive occurs at a point between the grid and the collector in Figure 4-2 (f). The space charge is more strongly felt on the right hand side of the grid, whereas on the left side, the potential merely curves upwards. In the region between points 1-4, two gaps in the Maxwellian distribution are introduced, again one due to the fact that there are no initial electron kinetic energies (at point 1) that are equal to 𝜓(𝑥). For this figure, we have 𝑣1 = (2 −𝜓(𝑥)𝑚)12and 𝑣2 = (2 𝜓𝑀−𝜓(𝑥)𝑚)12. The distribution function is equal to 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × { Θ(𝑣𝑥 − 𝑣1) − Θ(𝑣𝑥 + 𝑣1) + Θ(𝑣𝑥 + 𝑣2)}. For points 4-6, the electrons flow in both directions and the middle gap is absent. The velocity distribution function is given by 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 + 𝑣2). Electrons at points 6-8 already possess enough energy to surmount the barrier, 𝜓𝑀 and therefore move unidirectionally with velocity distribution function given by 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) ×Θ(𝑣𝑥 − 𝑣1). In Figure 4-2 (g), the space charge is strongly felt on the both sides of the grid, leading to two peaks. However, the peak in the region between the collector and the grid has a higher value, 𝜓𝑀 than the local extremum between the grid and the emitter, 𝜓𝑀,𝑙. In the region between points 1-3, the electrons move bidirectionally with the distribution function given by 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 + 𝑣1), where 𝑣1 = (2 𝜓𝑀−𝜓(𝑥)𝑚)12. At points 3-9, two areas are carved out of the Maxwellian distribution due to the foregoing reasons. The distribution function in this region is given by 𝑓(𝑥, 𝑣) =2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × { Θ(𝑣𝑥 − 𝑣1) − Θ(𝑣𝑥 + 𝑣1) + Θ(𝑣𝑥 + 𝑣2)}, where 𝑣2 =(2 𝜓𝑀,𝑙−𝜓(𝑥)𝑚)12. At points 9-11, the middle gap (– 𝑣2 ≤ 𝑣𝑥 ≤ 𝑣2 ) is removed and the distribution is the same as the one between points 1-3. Finally, at points 11-13, electrons can only move unidirectionally and 75 their distribution function is given by 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 −𝑣1). In Figure 4-2.(h), two motive peaks arise at the regions between the emitter-grid and the grid-collector. The peak situated between the grid and the emitter is at a higher potential than the local extremum point, 𝜓𝑀,𝑙, between the grid and the collector. At points 1-3, the electrons move bidirectionally with the distribution function 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 + 𝑣1), where 𝑣1 =(2 𝜓𝑀−𝜓(𝑥)𝑚)12. At points 3-13, electrons already possess enough kinetic energy to complete their journey towards the collector and therefore move unidirectionally with 𝑓(𝑥, 𝑣) =2 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑀−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) × Θ(𝑣𝑥 − 𝑣1). In the power generation mode of the TECs, the motive of the point outside the emitter will in most cases be below the motive of the point just outside the collector. This means that the statistics of the electrons is radically different in the presence of a grid and the conventional TEC equations can lead to physically untenable solutions. For instance, not carving out the middle part of the velocity distribution will make the concentration of the electrons in the negative motives (at the proximity of the grid) exponentially high; this upshot is an artifact of neglecting that electrons cannot have negative energies at the point just outside the emitter (Figure 4-2 (c, e, f, g)). Poisson equation The electron density, 𝑛(𝑥), is calculated by integration of the velocity distribution function as 𝑛(𝑥) =∫ ∫ ∫ 𝑓(𝑥, 𝑣) 𝑑𝑣𝑥 𝑑𝑣𝑦 𝑑𝑣𝑧∞−∞∞−∞∞−∞. The distribution function takes different forms depending on the motive profile, as fully detailed in the preceding section. The overall mathematical form of 𝑓(𝑥, 𝑣) can be captured by two different formulae. For the cases where 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑖) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32 exp (𝜓𝑖−𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸) ×Θ(𝑣𝑥 ± 𝑣1), where 𝑣1 = (2 𝜓𝑖−𝜓(𝑥)𝑚)12, and 𝑖 is the reference point in calculating the distribution functions, i.e., the point of maximum motive, the resulting electron density is 𝑛(𝑥) = 𝑛(𝑥𝑀) exp (𝜓𝑖 − 𝜓(𝑥)𝑘𝐵𝑇𝐸)(1 ± erf (√𝜓𝑖 − 𝜓(𝑥)𝑘𝐵𝑇𝐸)). (4-1) 76 Numerical difficulties can arise in calculation of the 1 − erf (√𝜓𝑖−𝜓(𝑥)𝑘𝐵𝑇𝐸), when the argument of the erf function is large—particularly at negative motive values—since the error function rapidly approaches 1 in these cases. We showed how this issue can be bypassed by using asymptotic expansion of the error function as in Chapter 3. The other possibility for the electron velocity distribution is given by 𝑓(𝑥, 𝑣) = 2 𝑛(𝑥𝑖) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32exp(𝜓𝑖 − 𝜓(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸)× { Θ(𝑣𝑥 − 𝑣1) − Θ(𝑣𝑥 + 𝑣1) + Θ(𝑣𝑥 + 𝑣2)} (4-2) , where 𝑣2 = (2 𝜓𝑗−𝜓(𝑥)𝑚)12 and 𝑣1 = (2 𝜓𝑖−𝜓(𝑥)𝑚)12 with 𝑥𝑗 as the position of the local extremum or the boundary point, provided the space charge does not create a peak in the inter-electrode region. The electron density in this case can be calculated by the integration of the velocity distribution function as 𝑛(𝑥) = 𝑛(𝑥𝑖) exp (𝜓𝑖−𝜓(𝑥)𝑘𝐵𝑇𝐸) (1 − 2 erf (√𝜓𝑗−𝜓(𝑥)𝑘𝐵𝑇𝐸) + erf (√𝜓𝑖−𝜓(𝑥)𝑘𝐵𝑇𝐸)). (4-3) Numerically calculating this equation can be particularly challenging for large values of the difference between the motive and the reference points, since the argument involves subtracting two numbers that are very close to 2 . We circumvented this issue by numerically integrating the velocity distribution function using the method of global adaptive quadrature [93]. The Poisson equation is solved in the two regions as the sum of the Laplace and Dirichlet Poisson equations with homogeneous boundary conditions as follows: 𝜓1 = 𝜓1,𝑑𝑝 + 𝜓1,𝑙 and 𝜓2 = 𝜓2,𝑑𝑝 + 𝜓2,𝑙, where subscripts 1 and 2 denote 0 < 𝑥 ≤ 𝑥𝐺 and 𝑥𝐺 < 𝑥 < 𝑑, respectively. 𝑥𝐺 corresponds to the position of the grid. The Laplace boundary condition is also set as 𝜓(𝑥𝐺) = 𝑒𝑉𝐺, where 𝑉𝐺 is the applied bias to the grid. Equation (4-3) is subsequently substituted into Poisson’s equation (1-6), subject to the boundary condition, 𝜓(0) = 0, 𝜓(𝑥𝐺) = 𝑒𝑉𝐺 and 𝜓(𝑑) = 𝑒 𝑉𝑖𝑛𝑡, where 𝑒 is the electron’s charge, 𝑉𝑖𝑛𝑡 is the voltage difference between the collector, and the emitter, 𝑑 is the inter-electrode distance and 𝜖0 is the permittivity of vacuum. 𝜓1,𝑑𝑝(𝑥) and 𝜓2,𝑑𝑝(𝑥) can be solved by discretization of the one-dimensional Laplacian 77 operator. Also, the solutions to the Laplace equation are 𝜓1,𝑙 = −𝑥𝑥𝐺(𝑒𝑉𝐺) and 𝜓2,𝑙 = −𝑒𝑉𝑖𝑛𝑡 +(𝑑 − 𝑥) (𝑒(𝑉𝑖𝑛𝑡−𝑉𝐺)𝑑−𝑥𝐺). Current density can be calculated as 𝐽 = 𝑒 ∫ ∫ ∫ 𝑣𝑥𝑓(𝑥, 𝑣) 𝑑𝑣𝑥 𝑑𝑣𝑦 𝑑𝑣𝑧∞−∞∞−∞∞−∞. As per continuity of charge, the value of the current density at steady state should be constant in the inter-electrode region. It is noted that in each one of the potential profiles in Figure 4-2 the value of the current density is the same at every point in the inter-electrode distance and equal to 𝑒 𝑛(𝑥𝑀) (2 𝑘𝐵 𝑇𝐸𝜋𝑚)1/2. The value of current density can be related to the saturation current density, 𝐽𝑠𝑎𝑡, obtainable from the Richardson-Dushman equation by observing that the current density can only reach its maximum value, i.e. 𝐽𝑠𝑎𝑡, when no extremum points occur in the inter-electrode region. In such cases, the value of the minimum velocity for the points just outside the emitter is equal to zero. In other words, the maximum current density is attained once the electron velocity distribution is HM at the point outside the emitter. Therefore, 𝐽𝑠𝑎𝑡 = 𝑒 ∫ ∫ ∫ 2𝑣𝑥 𝑛(𝑥𝑀) (𝑚2 𝜋𝑘𝐵 𝑇𝐸)32exp (𝜓𝑀𝑘𝐵 𝑇𝐸−𝑚𝑣22 𝑘𝐵𝑇𝐸)∞−∞∞−∞∞−∞× Θ(𝑣𝑥) 𝑑𝑣𝑥 𝑑𝑣𝑦 𝑑𝑣𝑧 = 𝐽 exp (𝜓𝑀𝑘𝐵 𝑇𝐸). (4-4) After solving the Poisson equation, the value of the maximum motive can be used to calculate the current density. Moreover, from solving the Vlasov equation, the value of current density can be related to the electron density at the point of maximum motive (i.e., 𝐽 = 𝑒 𝑛(𝑥𝑀) (2 𝑘𝐵 𝑇𝐸𝜋𝑚)1/2). The self-consistent strategy to solve the Poisson and Vlasov equations is presented in the flowchart in Figure 3-1. The electron density at the point of maximum motive is calculated by assuming that the current density is equal to the saturation current density in the first approximation. Based on these and considering the different possibilities that arise in Figure 4-2, the electron velocity distribution at different regions are calculated. These functions are subsequently integrated to calculate the electron density. The electron density is substituted inside the Poisson equation to calculate the potential profile. The value of the current density for the next iteration is calculated from equation (4-4) and the solution of the Poisson equation (to determine the point of maximum motive). At the end of each iteration, it is checked whether the value of the calculated current density is in accordance with the current density from the previous step. 78 Results and discussion The self-consistent approach, as depicted in Figure 3-1, was used to calculate the motive and the electron density in the inter-electrode region of TECs in the presence of a grid at different locations. The proposed model was implemented in MATLAB. At each iteration, after solving the Poisson equation and updating the motive, the program checks for proper statistics that applies to the electrons. It is possible to encounter different types of statistics at different iterations, since the position of the extremum points vary depending on the electron density. Moreover, in updating the motive at each iteration, the value of the new potential calculated from the Poisson equation is mixed with the motive from the previous step, in order to avoid the strong fluctuations that arise in solving this highly non-linear system of equations. This is formulated as 𝜓𝑛(𝑥) = 𝛼𝜓ℎ𝑝(𝑥) + (1 − 𝛼)𝜓𝑛−1(𝑥), where 𝜓𝑛(𝑥) represents the updated motive before addition of the Laplace solution, 𝜓𝑛−1(𝑥) represents the motive before addition of the Laplace solution from the previous iteration, and 𝜓ℎ𝑝(𝑥) is the solution to the Poisson equation with the Dirichlet boundary conditions. 𝛼 is a number between 0 and 1. For cases where convergence can be easily attained, a value of 𝛼 = 0.01 would be adequate. In circumstances where wild fluctuations in the motive occur, 𝛼 = 0.001 is more apt. It is assumed that the grid’s length scale is finer than the Debye length of the electron gas (typically 10 − 500 𝜇𝑚 in TECs). Under this condition, the plane of the grid (including the holes) can be considered to be an equipotential surface and thus the motive at the gate position can be set equal to 𝑒 𝑉𝐺, where 𝑉𝐺 is the grid bias with respect to the emitter. The Poisson equation is solved recursively with the Vlasov equation using the strategy portrayed in Chapter 3. At each iteration, the solution is mixed with the solution from the previous step to avoid oscillations. The mixing ratio between the newly acquired motive in the current iteration and the previous motive, 𝛼, is chosen a small number, e.g., 0.001, in the beginning since the solutions to the Vlasov and Poison equations are highly decoupled. This value is subsequently increased to 0.01 to expedite the convergence of the solutions (Figure 4-3.(a)). The convergence criterion is chosen such that the cumulative root mean square change in motive is less than 0.1 % .Figure 4-3.(a-b) depict the evolution of the motive and the electron density as a function of the iteration number. The coupled Poisson-Vlasov system is such that an increase in one variable results in its reduction in the subsequent step. For instance, if the potential increases, the electron density is significantly reduced from the Vlasov equation, leading to a smaller potential in the following iteration. The presence of velocity gaps, mentioned above, are extremely important in this regard. If the solutions to the Vlasov equation to a traditional TEC are employed (without the gaps in the phase-space), the electron density significantly increases in negative motives, making a simultaneous solution to both equations unattainable. This can be seen in this fashion: if 79 the regular TEC equations are employed, the electron density in the potential well is determined by the 𝑛(𝑥) = 𝑛(𝑥𝑀) exp(𝛾) (1 + erf 𝛾1/2). This means that at negative motives, the electron density grows exponentially higher than 𝑛(𝑥𝑀,𝑙) with distance. This reduced electron density at 𝑥𝑀,𝑙 increases the concentration of electrons that can overcome the barrier in the next iteration. Therefore the subsequent electron density in the next iteration will be substantially higher and this goes on ad infinitum making solutions to the system unattainable. However, as depicted in Figure 4-3. (c), the electron density in the region 𝑥 ∈ (0, 𝑥𝐺) is a decreasing function of distance which guarantees the convergence of the solutions. In a general case, a more rigorous argument can be obtained by considering the electron density equation in this region. The maximum value of 𝛾𝑙 for which 2 erf 𝛾𝑙1/2≤ erf 𝛾1/2 is equal to 0.4769. Therefore, whenever 𝛾𝑙 plunges to a value smaller than this, 𝑛(𝑥) would be an increasing function of distance in that region. However, the convergence can be rescued by considering that in the next iteration, the value of 𝜓𝑀,𝑙 increases (since more electrons are now allowed to pass the local barrier), and thereby the value of 𝛾𝑙 enters the converging region. Any two coupled partial differential equations that are compatible (in the sense that a change in one does not amplify in the other) can be such that a common solution two both of them does not exist. In such cases, one would expect an oscillatory response. The lack of oscillatory behaviour in our cases was checked by setting the mixing ratio to a higher value (𝛼 = 0.1) after the convergence was reached (𝑛 = 1000 in Figure 4-3. (a)) and carrying out the calculations for another 1000 iterations. The solutions remained constant in all the cases that were studied in this Chapter, denoting a non-oscillatory response. 80 Figure 4-3. (a) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.001 for 𝑛 ≤ 500 and 0.01 for 500 < 𝑛 < 1000 and 0.1 for 𝑛 > 1000. The workfunctions are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, the emitter’s temperature, 𝑇𝐸, is 1800 𝐾. The collector’s temperature is assumed to be low enough so that the back-emission can be neglected. The inter-electrode gap distance, 𝑑, is 100 𝜇𝑚 and the grid is situated at 𝑥𝐺 = 10 𝜇. The bias on the grid, 𝑉𝐺 , is 2 V. The applied voltage between the collector and emitter is 0 𝑉. b)The convergence of the motive in the IER for 81 different values of the iteration number, 𝑛. (c) Evolution of the electron density as a function of position for different iterations. Figure 4-4. The motive profile and the electron density of a TEC with 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, 𝑇𝐸 = 2,200 𝐾, 𝑑 = 100 𝜇𝑚 , 𝑥𝐺 = 10 𝜇𝑚 and 𝑉 = 0 𝑉. Also shown are the potential profile and the electron density for the same TEC, but in the absence of a grid. As detailed in the main text, the overall effect is a total reduction in the maximum motive barrier. The motive and electron density distributions of a TEC are depicted in Figure 4-4. The parameters of the TEC are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, 𝑇𝐸 = 2,200 𝐾, 𝑑 = 100 𝜇𝑚 , 𝑥𝐺 = 10 𝜇𝑚 and 𝑉 = 0 𝑉. The position of the grid is chosen such that the tallest motive peak occurs after the grid, in accordance with the conditions necessary to generate the momentum gaps. The presence of velocity gaps, mentioned above, are indispensable to the stable operation of this device. In their absence, the system becomes unstable and therefore not physically realizable. This can be understood by the following argument. For this transport system to be physically possible, there needs to be a negative feedback mechanism between the concentration of the electrons and the motive profile that 82 is generated due to the presence of these electrons and the external potentials. This requirement is rooted in the fact that electrons are only ejected from the emitter, and henceforth their concentration needs to decrease in the region between point 2 and the grid (Figure 4-1), otherwise they are being accumulated inside the well. However, this cannot occur in steady state mode of operation of the device. This argument is supported mathematically by noting that the electron density function in the potential well region in Figure 4-2 is a decreasing function of distance, except for when the value of 𝛾𝑙 plunges below 0.4769. However this does not lead to instability since the local barrier height on the left will increase due to accumulation of charges in the well, which would in turn increase the value of 𝛾𝑙 , rendering the system stable. Figure 4-5. (a) The converged motive profile of a TEC with 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, 𝑇𝐸 = 3,000 𝐾, 𝑑 = 1 𝑚𝑚, 𝑥𝐺 = 500 𝜇𝑚 for different values of the grid bias. (b) the resulting electron density as a function of position for different gate biases. (c) the maximum motive, 𝜓𝑀 as a function of the grid bias. The computed data can be estimated by a double exponential fit, 𝜓𝑀 = 𝑎 𝑒𝑥𝑝(𝑏 𝑉𝐺) + 𝑐 𝑒𝑥𝑝(𝑑 𝑉𝐺), where 𝑎 = 0.03798, 𝑏 =−3.983 , 𝑐 = 1.075 and 𝑑 = −0.2363 with 𝑅𝑀𝑆𝐸 = 0.01705. 83 The motive and the electron density vs. distance as a function of the grid bias are depicted in Figure 4-5. (a-b) for a TEC with the parameters given in the figure caption. The maximum motive as a function of the applied grid voltage is plotted in Figure 4-5.(c). It is interesting to note that increasing the grid bias impacts the device performance in two opposing ways: the overall electron density is increased by incrementing the grid bias according to Figure 4-4 and Figure 4-5.(b), since more electrons are allowed to be in the inter-electrode region compared to the case that no grid is present. This is due to the presence of higher potential energy differences, leading to a larger spectrum of kinetic energies and therefore velocities. Therefore, the higher number of electrons in the inter-electrode region increases the space charge potential. The other effect is that the motive boundaries (the electrodes and the grid) are shifted downwards due to the negative applied bias. Mathematically, the first effect is due to the Vlasov equation, whereas the second one is due to the Poisson equation. The interplay between these factors leads to an overall reduction in the maximum motive (Figure 4-4 and Figure 4-5.(c)) and hence an increase in the current density. Figure 4-6. The maximum motive (red curves) and the current density (blue curves) as a function of the grid bias for different values of 𝑑. The current densities in the absence of the grid are plotted as straight lines. The parameters of the TEC are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, 𝑇𝐸 = 1,800 𝐾 and 𝑥𝐺 = 𝑑/2. The outcome is that the solution to the Poisson equation has a smaller peak compared with the case that no grid is present (Figure 4-4 and Figure 4-5. (a), (c)). The presence of these two effects can be demonstrated in Figure 4-5.(c), where it is noted that a single exponential curve cannot capture the trend of 84 changes of 𝜓𝑀, and that a bi-exponential function is necessary for this purpose. It is also noted from Figure 4-5. (c) that even setting 𝑉𝐺 = 0 is more favorable than not having a gate at all. The developed model can be used to calculate the value of the grid bias necessary to obtain a certain amount of current density or motive barrier. The results are depicted in Figure 4-6 for a TEC (see the device parameters in the figure caption) at different values of 𝑑. It is observed that at smaller values of 𝑑—mostly applicable to emergent TECs—smaller values of grid bias, on the order of several volts can diminish the space charge effect and improve the output power density by an order of magnitude. At higher distances, more bias needs to be applied to the grid according the Figure 4-6 to obliterate the space charge, since more electrons inhabit the inter-electrode region. Remarkably, the current densities in all of the three cases studies in Figure 4-6 are smaller when no grid voltage is applied. Figure 4-7.(b) demonstrates the convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛. This TEC is comprised of two electrodes with the same workfunction of 𝜙𝐸 =𝜙𝐶 = 4.5 𝑒𝑉, the emitter’s temperature, 𝑇𝐸, is 3,000 𝐾—this high temperature is chosen deliberately to generate a substantial amount of space charge so that different motive schemes can be easily observed. The collector’s temperature is assumed to be low enough so that back-emission of current from the collector is negligible. The inter-electrode gap distance, 𝑑, is 1 𝑚𝑚 and the grid is situated half-way between the emitter and the collector, at 𝑥𝐺 = 500 𝜇𝑚. A positive bias of 2 𝑉 is applied to the grid with respect to the emitter (𝑉𝐺). 85 Figure 4-7 (a) The convergence of the motive in the inter-electrode region for different values of the iteration number, 𝑛. The workfunctions are 𝜙𝐸 = 𝜙𝐶 = 4.5 𝑒𝑉, the emitter’s temperature, 𝑇𝐸, is 3,000 𝐾. The collector’s temperature is assumed to be low enough so that the back-emission can be neglected. The inter-electrode gap distance, 𝑑, is 1 𝑚𝑚 and the grid is situated half-way between the emitter and the collector, at 𝑥𝐺 = 500 𝜇. The bias on the grid, 𝑉𝐺 is 2 V. The applied bias between the collector and emitter is −2 𝑉. (b) Evolution of the electron density as a function of position for different iterations. c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.0001 for 𝑛 ≤ 8000 and 0.001 for 1200 > 𝑛 > 800 and 0.01 for 𝑛 > 1200. The internal voltage, which is equal to the applied bias between the collector and the emitter in this case, is −2 𝑉. In the first iteration, the potential is linear since the initial concentration of electrons is assumed to be zero (linear motive). Therefore, the point of maximum motive coincides with the point outside the collector. Here, the point with 𝜓(𝑥) = 0 is the point right outside the emitter. The electron density in the inter-electrode region is plotted in Figure 4-7.(b) for different values of iteration number. In the beginning, since no peaks are present in the inter-electrode region, all electrons that are emitted from the emitter can be present in the inter-electrode region. Therefore, for initial values of 𝑛, 86 the electron density is high. However, in later iterations, this high number of electrons produces peaks in the inter-electrode region. These peaks, reduce the number of electrons that have enough kinetic energy to overcome them, and therefore their concentration is reduced in the forthcoming iterations. Eventually an equilibrium will be reached where the motive and the velocity distributions are simultaneous solutions of the Poisson and Vlasov equation. The value of current density as a function of the iteration number is presented in Figure 4-7.(a). During the first iterations, the current density is very close to its maximum value, i.e., the saturation current density, since no peaks are present. As 𝑛 increases, the presence of the peak in the motive reduces the current density. This reduced current density in turn leads to a lower number of electrons and therefore a smaller motive peak. These opposing forces balance at higher values of 𝑛, when equilibrium is reached. When the point outside the collector and the emitter lied at the same motive, the electron statistics of the device transitioned from Figure 4-2.(g) to Figure 4-2.(h). This behaviour is illustrated in Figure 4-8.(a). The parameters of this TEC are the same as the TEC in Figure 4-7. The evolution of the motive as a function of distance for different iterations is illustrated in Figure 4-8.(a) The same interplay between the motive and electron density is observed. Notably, the high value of current density (or electron density) at early stages, leads to high values of the motive peak, which subsequently reduce the electron density. As depicted in Figure 4-8.(b), the motive reaches a high peak at 𝑛 = 100, which is subsequently relaxed to the equilibrium value of the potential at higher values of 𝑛. The evolution of the electron density, 𝑛(𝑥), as a function of the distance for different iterations is plotted in Figure 4-8.(c). As expected, the values of the electron density start up at a high value and reach equilibrium after fluctuations. The evolution of the current density as a function of the potential is presented in Figure 4-8.(a) 87 Figure 4-8. (a) The convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛, for a TEC with similar parameters as that of Figure 4-7. The applied bias between the collector and emitter is 0 𝑉. (b) Evolution of the electron density as a function of position for different iterations. (c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.0001 for 𝑛 ≤ 8000 and 0.001 for 1200 > 𝑛 > 800 and 0.01 for 𝑛 > 1200. When the internal motive between the emitter and the collector is negative, i.e. when a positive bias is applied between the emitter and the collector, assuming that the workfunctions of the two electrodes are equal, the initial current density at 𝑛 = 1 will be equal to the saturation current density. This is in contrast to Figure 4-7, when the initial value of the current density was equal to 𝐽 = 𝐽𝑠𝑎𝑡 exp (−|𝑒 𝑉𝑖𝑛𝑡|𝑘𝐵 𝑇𝐸). The motive profile as a function of distance for different iteration numbers is presented in Figure 4-9.(b) This motive belongs to a TEC with a grid in its inter-electrode region with similar parameters as that of Figure 4-7. The internal bias in here is 𝑉𝑖𝑛𝑡 = 2 𝑉. The electron density as a function of position for different iterations and the current density as a function of the iteration number are presented in Figure 4-9.(c)and Figure 4-9.(a), respectively. 88 Figure 4-9. (a) The convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛, for a TEC with similar parameters as that of Figure 4-7. The applied bias between the collector and emitter is 2 𝑉. (b) Evolution of the electron density as a function of position for different iterations. (c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.0001 for 𝑛 ≤ 8000 and 0.001 for 1200 > 𝑛 > 800 and 0.01 for 𝑛 > 1200. The motive profile is strongly dependent on the position of the grid (Figure 4-10). For instance, when the grid is placed half-way between the emitter and the collector, the maximum motive occurs at the grid-collector side. However, when the grid is placed at 𝑥𝐺 =910𝑑, the peak will be shifted towards the emitter-grid side. It is noted that before reaching an equilibrium point, the motive spans a wide array of the possible scenarios of the motive profiles given in Fig 2. In the beginning (𝑛 = 1), the motive profile is similar to that of Figure 4-2.(c) As the simulation progresses, a peak on the emitter-grid side is generated and the statistics become similar to that of Figure 4-2.(e). Subsequently, at later iterations, according to Figure 4-10.(b), this peak’s height surpasses the maximum motive’s height on the grid-collector (i.e., the point right outside the collector). Therefore, the electron statistics of this TEC transition from Figure 4-2.(e) to Figure 4-2.(h). This peak is reduced later on and the electron statistics is reverted to Figure 4-2.(e). The 89 evolution of the electron density and the current are illustrated in Figure 4-10.(b) and Figure 4-10.(c), respectively. Figure 4-10. (a) The convergence of motive in the inter-electrode region for different values of the iteration number, 𝑛, for a TEC with similar parameters as that of Figure 4-9. The applied bias between the collector and emitter is −2 𝑉. The grid is positions in 𝑥𝐺 =910𝑑 and 𝑉𝐺 = 10 𝑉 (b) Evolution of the electron density as a function of position for different iterations. (c) Evolution of the emitter’s current density, 𝐽𝐸, as a function of the iteration number. The value of 𝛼 is chosen as 0.001 for 𝑛 ≤ 1,000 and 0.01 for 𝑛 > 1000. The grid loss can be easily included in this model by multiplying the current density at each iteration level by (1 − 𝑙), where 𝑙, represents the loss fraction of the grid. This reduced current density leads to a lower 𝑛(𝑥𝑀) and therefore lower space charge at the expense of power loss through the grid. For all the cases studied here, the grids are assumed to be lossless. This model neglects the electronic structure of the grid metal and its impact on the energy dependence of the absorption/reflection/transmission coefficients of the impinging electrons; for cases where this might become important, such as in atomically thin grids, the model would need to be modified to incorporate this energy-dependent behaviour. 90 The work due to Meir et al. presents that by applying a longitudinal magnetic field, the transparency of a grid to electrons can be enhanced to their geometrical transparency (98 % in their case) [8]. However, by using nanomaterials, it is possible to improve the transparency even further [73]. If a bias of 2 V is applied to the grid, this geometric loss corresponds to ~ 4 % loss in the total output power if the space-charge is completely obliterated. On the other hand, if higher voltages are applied, it is necessary to reduce the grid loss further to improve the overall efficiency. Summary In summary, it was shown that momentum gaps could arise in the phase-space of electrons undergoing transport under certain conditions in vacuum. The transport model for this phenomenon was developed and was employed to analyze the performance of thermionic energy converters in the presence of grids. It was further demonstrated that the application of positive voltages to the grid can significantly mitigate the space charge effect. A model was proposed to analyze the output characteristics of a thermionic energy converter in the presence of a grid. It was shown that the electron statistics can be significantly different in the presence of a grid compared with a simple thermionic converter in the absence of the grid. The root of these differences were attributed to the fact that electrons originate from the emitter, and that the minimum initial kinetic energy of the electrons is zero. These modifications were considered in solving the Vlasov equation. The resulting velocity distribution function was used to calculate the electron density that was successively used to calculate the potential through the Poison equation. The Poisson-Vlasov system was solved self-consistently until divergence. This method was subsequently employed to calculate the output characteristics of thermionic converters in the presence of grids under different operating conditions. 91 Chapter 5- Low-pressure plasma-enhanced behaviour of thermionic converters5 The widely accepted method of reducing space-charge in thermionic converters is by introduction of positively-charged ions to neutralize the negative space charge of electrons. Although the behaviour of thermionic converters including high pressure plasmas has been studied extensively, the physics of their operation under low pressure conditions is less well understood. However, due to current advances in nanotechnology such as the possibility to intercalate nanomaterials-based electrodes with alkali metals to reduce their workfunction with lower pressures of the neutral alkali metals in the inter-electrode region, low-pressure plasma energy converters may provide a solution to some of the long-standing challenges such as the space charge effect. Here, we develop the physics of low-pressure thermionic converters by providing closed-form solutions to a Vlasov-Poisson system. We demonstrate that a host of possibilities could arise due to the intricate interaction between spatially varying electron and ion concentrations, leading to phenomena such as plasma oscillations at higher ion fluxes. We show that relatively low flux densities of ions (~ 5 × 10−4 that of the flux density of electrons) could improve the electron current density by a factor of 7. We further extend the model by including electron and ion emission from both the cathode and the anode electrodes. Theory and model Depending on the concentrations of the electrons and ions, a whole variety of different motive profiles could arise in the inter-electrode region. The motives differ in the number of peaks and valleys and the magnitude of each peak with respect to other peaks. It is assumed that ions have a unit positive charge and therefore the motive that they experience is the negative of the electron motive. The motives start at the point outside the emitter (𝑥 = 0) and end at the point outside the collector (𝑥 = 𝑑). Two special cases of motive profiles with 3 valleys and 2 peaks in the inter-electrode region are depicted in Figure 5-1. These motive profiles can arise both for electrons and ions (with the motive for positive ions being the negative of that for electrons) depending on their concentrations. 5 This chapter is reproduced with minor modifications from [5]. 92 We will describe the species velocity distribution in these cases and then provide a general algorithm that can be applied to an arbitrary motive profile. For points 3-4 and 11-12 in the motive on the left, the velocity distribution functions are the same as that of regular TECs, and have been described before [2]. The velocity distribution functions are 𝑓𝛼(𝑥, 𝑣) = 2 𝑛𝛼(𝑥𝑀) (𝑚𝛼2 𝜋𝑘𝐵 𝑇𝐸)32exp (𝜓𝑀,𝛼 − 𝜓𝛼(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝛼𝑣22 𝑘𝐵𝑇𝐸)× Θ(𝑣𝑥 ± 𝑣𝑥,𝑚𝑖𝑛), (5-1) where 𝛼 ∈ {𝑒, 𝑖}, with 𝑒 and 𝑖 representing electrons and ions, respectively. 𝑛𝛼(𝑥𝑀) is the species 𝛼’s concentration at the point of maximum motive 𝜓𝑀,𝛼. 𝑚𝛼 is the mass of species 𝛼, 𝑇𝐸 is the emitter’s temperature, 𝑘𝐵 is the Boltzmann constant, 𝑥 is the direction of propagation, 𝑣 is the speed and 𝑣𝑥 is the velocity along the 𝑥 direction. Also, 𝜓𝑀,𝑖 = −𝜓𝑀,𝑒, and 𝑣𝑥,𝑚𝑖𝑛 represents the minimum permissible velocity of the charge carrier at position 𝑥, and is equal to 𝑣𝑥,𝑚𝑖𝑛 = (2 𝜓𝑀,𝛼−𝜓𝛼(𝑥)𝑚𝛼)12. Θ is the unit step function, and the positive sign in its argument is for 𝑥 ≤ 𝑥𝑀 (since negative velocities are permissible) and the negative sign for 𝑥 > 𝑥𝑀 [14]. In order to obtain these equations, the 1-D Vlasov equation is solved assuming that the electrons and ions adopt a hemi-Maxwellian distribution at the point of their maximum motive (different points for electrons and ions). Equation (5-1) is similar to equations (2-1a) and (2-1b), but is more general in the sense that it can applied to both electrons and ions. However, when the charge carriers are in regions 1-3 and 4-11, a gap arises in their velocity distribution function due to the reflection of the charge carriers from the local maxima: 𝑓𝛼(𝑥, 𝑣) = 2 𝑛𝛼(𝑥𝑀) (𝑚𝛼2 𝜋𝑘𝐵 𝑇𝐸)32exp (𝜓𝑀,𝛼 − 𝜓𝛼(𝑥)𝑘𝐵 𝑇𝐸−𝑚𝛼𝑣22 𝑘𝐵𝑇𝐸)× { Θ(𝑣𝑥 − 𝑣1) − Θ(𝑣𝑥 + 𝑣1) + Θ(𝑣𝑥 + 𝑣2)}, (5-2) where 𝑣2 = (2 𝜓𝑀,𝑙,𝛼−𝜓𝛼(𝑥)𝑚𝛼)12 and 𝑣1 = (2 𝜓𝑀,𝛼−𝜓𝛼(𝑥)𝑚𝛼)12. 𝜓𝑀,𝑙,𝛼 is the local minimum (including the boundary points), and in region 1-3 is the motive of the emitter, i.e., 0, whereas in the region 4-11 it is equal to 𝜓4. Therefore, for the velocity distributions at the top of Figure 5-1, the local minima of motive to calculate velocities 𝑣2 and 𝑣3 are given by the motive at points 1 and 4, respectively. These distribution functions can be obtained based on the arguments that we made in [submitted] to analyze the electron dynamics in the presence of a grid. Briefly put, since charge carriers are assumed to originate from one electrode and travel to the other, in each well, only electrons with initial kinetic energies higher than all the preceding peaks, since otherwise they would not be able to reach that region and would be reflected back. 93 In contrast to the situation with a grid, in the case of the low-pressure plasma device, momentum gaps could arise in multiple regions of space. For instance, for the motive on the right hand side of Figure 5-1, 3 regions with distinct momentum gaps arise. These regions are distinguished by having different values of maximum local motives. Figure 5-1. The motive profile and the resulting distribution functions for electrons in two special cases where the number of valleys and peaks are less than 3 (the motive profile for positive ions is the negative of that of electrons). 𝜓 is the motive; 𝑓(𝑥, 𝑣𝑥) is the velocity distribution function. 𝑣𝑥 represents the velocity along the direction of propagation (𝑥). Point 1 corresponds to the point just outside the emitter and the last point (12 or 11) represents the point just outside the collector. In the left motive, 3 regions with different velocity distribution functions exist, and there are two regions that have two different momentum gaps. For the motive on the right, there are 4 distinct regions and 3 different momentum gaps. These motives are due to the combined effects of both electrons and ions. For the velocity distributions at the top, the local minima of motive to obtain velocities 𝑣2 and 𝑣3 are given by the motive at points 1 and 4, respectively, whereas for the velocity distributions at the bottom, the local minima of motive to obtain velocities 𝑣2, 𝑣3 and 𝑣4 are given by the motive at points 1, 4 and 7, respectively. For the velocity distributions at the bottom of Figure 5-1, the local minima of motive to obtain velocities 𝑣2, 𝑣3 and 𝑣4 are given by the motive at points 1, 4 and 7, respectively. Although in Figure 5-1, the number of the reported peaks and valleys are below 4, equations (5-1) and (2) can be modified to be applicable to any number of peaks and valleys. Assume that the peaks in the motive are numbered from 1 to 𝑛𝑝, where 𝑛𝑝 is the total number of peaks. At an arbitrary point 𝑥, the 94 distribution function is determined by the positive sign in the argument of equation (5-1) (electrons traveling unidirectionally), if 𝑥 lies on the right hand side of peak 𝑛𝑀, where 𝑛𝑀 is the index for the peak with the highest motive. The distribution function is determined by the negative sign in the argument of equation (5-1) for all points lying to the left hand side of the first peak, 𝑛1, or at points between the peaks 𝑛𝑖 and 𝑛𝑗 (𝑗 > 𝑖) where the motive is higher than the motive of peak number 𝑛𝑖 and lower than the motive at peak number 𝑛𝑗. At all other points x, the distribution is determined by equation (5-2) with different values of 𝑛𝛼(𝑥𝑀), 𝜓𝑀,𝑙,𝛼 and 𝜓𝑀,𝛼 determined by the highest peak prior to (to the left of) point 𝑥. The resulting electron distributions are integrated and inserted into Poisson’s equation, 𝑑2𝜓(𝑥)𝑑𝑥2=−𝑒2𝜖0(𝑛𝑒(𝑥) − 𝑛𝑖(𝑥)), where 𝑛𝛼(𝑥) represents the concentration of species 𝛼. The Poisson and Vlasov equation are solved recursively based on the method that we have proposed before [3]. Results and discussion When the flux of the ions is sufficiently low, their influence on the behaviour of the TEC is the total reduction of space-charge, with no dramatic changes in the motive profile. For instance, when the saturation flux density of the ions, 𝐽𝑠𝑎𝑡,𝑖, is 10−4 times smaller than the saturation current density of the electrons, 𝐽𝑠𝑎𝑡,𝑒, (given by the Richardson-Dushman equation [3]), i.e., when 𝛽 ≡𝐽𝑠𝑎𝑡,𝑖𝐽𝑠𝑎𝑡,𝑒= 10−4, the maximum motive is reduced by 0.05 𝑒𝑉 (Figure 5-2 .(a)) in a TEC with the parameters given in the caption of Figure 5-2. For each of the cases studied in Figure 5-2, a mixing ratio parameter, 𝛿, ([3]) is defined which determines the weight by which the motive is updated in the carrier density at each iteration. Each iteration yields the current density obtained from the solution to the Vlasov equation, i.e., the steady state current density obtained if the motive is frozen at its given value. The motive and carrier densities in Figure 5-2 all converge and, therefore, it can be concluded that the system is stable and reaches a DC steady-state (no plasma oscillations arise). Figure 5-2. (c) shows the electron density as a function of distance. It can be observed that in the presence of ions, a higher number of electrons are allowed in the inter-electrode region, although the overall motive is reduced in Figure 5-2. (a). Moreover, the electron density becomes smaller as the electrons propagate towards the collector, as one would expect in regular TECs as well. The ion density, on the other hand, has a minimum value somewhere in the inter-electrode region and has its highest values closer to the electrodes (Figure 5-2. (c)) when 𝛽 = 10−4. The situation can become different when 𝛽 = 5 × 10−4. In this case, the concentration of ions (Figure 5-2. (c)) is high enough to cause negative motives at some points (Figure 5-2. (a)). The concentration of electrons is nonetheless a decreasing function of distance, since they face only one barrier on their way to 95 the collector (Figure 5-2 .(a)). Moreover, since the current density plateaus at higher values of the iteration, the system is stable and a DC solution to this plasma enhanced TEC can be obtained. According to Figure 5-2. (a), the maximum motive is 0.24 𝑒𝑉 reduced compared to the case that no ions are present, corresponding to a 7-fold enhancement in output current density. If the saturation flux density of the ion emission is increased further, plasma oscillations could occur and the system will not have a stable DC behaviour in steady state. This phenomenon could only arise if the maxima of the electron and ion densities both occur at some point in the inter-electrode region, rather than the points just outside the emitter and collector. In such cases, the electrons tend to travel from the point of their highest concentration to a location with the highest concentration of ions. Therefore, a static steady state can never be attained and the electron and ion densities, as well as the motive, change cyclically and in accordance with each other. The only conceivable exception is when the point with maximum electron density coincides with the point with maximum ion density, which cannot occur based on the statistics described by equations 5-(1:2). In terms of motive, this phenomenon could arise when there are momentum gaps present in the phase-space of the propagation of both the electrons and the ions, so that there is at least one region for each species where the velocity distribution function follows equation (5-2). In the TEC described in Figure 5-2, when 𝛽 = 10−4, no momentum gap is present in the phase space of either electrons or ions. In contrast, when 𝛽 = 5 × 10−4, momentum gaps exist for the transport of ions, but not electrons. 96 97 Figure 5-2. The characteristics of a TEC in the presence of ions 𝛽 = 0, 10−4 and 5 × 10−4 . The workfunctions of the emitter and the collector are 𝜙𝐸 = 𝜙𝐶 = 2.5 𝑒𝑉, the emitter’s temperature is 𝑇𝐸 =1500 𝐾 and the inter-electrode distance is 𝑑 = 100 𝜇𝑚. The resulting (a) motive,𝜓(𝑥) (b) electron density, 𝑛𝑒(𝑥), and (c) ion density, 𝑛𝑖(𝑥), are plotted. The case with 𝛽 = 0 represents the absence of ions. Figure 5-3 portrays the characteristics of a TEC with a higher amount of ion saturation flux density (𝛽 = 10−3 of the electron saturation current density). Figure 5-3 (a) depicts the evolution of the current density as a function of the iteration number. The values of 𝛼 are typically chosen to be small in the beginning (𝛼 = 10−3) and then gradually increased to 𝛼 = 0.1 when the solutions to the Poisson and Vlasov equation become more consistent as the initial numerical instabilities have abated. The values are 𝛼 are reduced again to 𝛼 = 10−3 at higher stages. As observed from Figure 5-3. (a), a DC steady state current density cannot be obtained as the current oscillates between several values. According to Figure 5-3. (c) and Figure 5-3. (d), both electrons and ions have their maximum density in the inter-electrode region. The motive profile at 𝑛 = 2763 corresponds to the case where the electrons do not possess a momentum gap in their phase space, whereas the ions, seeing the negative of the motive experienced by the electrons, do indeed possess a momentum gap and will have their maximum concentration in the inter-electrode region. At 𝑛 = 3434, a new peak emerges and, beyond 98 this iteration, the electrons will adopt a momentum gap. 𝑛 = 3878 is the mid-point between the two more extreme cases of 𝑛 = 1294 and 𝑛 = 1300. This process repeats itself, leading to oscillating species densities and motives. Based on Figure 5-3. (a), a common solution to the Poisson and time-independent Vlasov equations does not exist, therefore, a DC steady state cannot be attained for this value of ions flux density. The only possibility is that a time-independent solution does not exist and therefore the system will have plasma oscillations. A time-dependent solution is beyond the scope of this work. Moreover, such cases are extremely rare in low-pressure plasma TECs. For instance, to obtain the ion flux in Figure 5-3, by using the potassium vapor pressure data [77] and the Saha-Langmuir equation [13], a potassium reservoir temperature of 870 𝐾 is required. This corresponds to a neutral potassium pressure of ~ 13 𝑘𝑃𝑎. The mean free path for the electrons and the gas molecules under such conditions can be approximated as ~ 70 𝜇𝑚 which is already smaller than the inter-electrode distance (100 𝜇𝑚), which is detrimental to device operation. Therefore, such cases are rarely subsumed within the operational regime of low-pressure plasmas. Nonetheless, an estimate for one of the time scales involved in the system may be obtained by noting that the velocity of the ions is on average several hundred times smaller than the velocity of electrons. For instance, the average velocity of potassium ions is ~ 270 times smaller than that of electrons. Therefore, this time scale is determined by the slower moving ingredients, i.e., ions. This can be estimated to be ~ 2𝑑 (𝑚𝛼𝑒 𝜓𝑀)1/2= 700 𝑝𝑠. This time scale is 2 orders of magnitude larger than the natural plasma oscillation period that can be found using the ion density. 99 100 Figure 5-3. The characteristics of a TEC with similar parameters as that of Figure 5-2 but in the presence of ions with the flux density of ions being 1 × 10−3 of that of electrons. These parameters lead to oscillations in the output characteristics of this TEC. Evolutions (in iteration number) of the (a) current density from the emitter, 𝐽𝐸, (b) motive, 𝜓(𝑥) (c) electron density, 𝑛𝑒(𝑥), and (d) ion density, 𝑛𝑖(𝑥), are plotted. The values of 𝛼 are depicted on the right axis in (a). Our model can be extended to include more general cases where back-emission of the ions and electrons from the collector is not negligible. Figure 5-4 shows such a case for a TEC with the parameters that are provided in the figure caption. The same strategy as before was applied to calculate the output characteristics. The charge density inserted into the Poisson equation was calculated as the sum of all contributors, i.e., the electron density from the emitter, 𝑛𝐸,𝑒(𝑥), ion density from the emitter, 𝑛𝐸,𝑖(𝑥), electron density from the collector, 𝑛𝐶,𝑒(𝑥), and ion density from the collector, 𝑛𝐶,𝑖(𝑥). Each one of these charge densities was calculated by integrating the corresponding velocity distribution function over the possible ranges of velocities. 101 102 Figure 5-4. The motive and carrier densities of a TEC with 𝜙𝐸 = 2.5 𝑒𝑉, 𝜙𝐶 = 2.0 𝑒𝑉, 𝑇𝐸 = 1500 𝐾, 𝑇𝐶 =1300 𝐾, and 𝑑 = 100 𝜇𝑚. The ratio of the ion to electron flux density is set equal to 10−4 in both collector and emitter. The (a) motive and the total species and the total absolute density 𝑛𝑡𝑜𝑡 ≡|𝑛𝐸,𝑒(𝑥) + 𝑛𝐶,𝑒(𝑥) − 𝑛𝐸,𝑖(𝑥) − 𝑛𝐶,𝑖(𝑥)|, (b) electron density resulting from the emitter, 𝑛𝐸,𝑒(𝑥) and collector, 𝑛𝐶,𝑒(𝑥) and (c) ion density resulting from the emitter, 𝑛𝐸,𝑖(𝑥) and collector, 𝑛𝐶,𝑖(𝑥), as a function of distance are plotted. The proposed method can be used to calculate the output current density-voltage characteristics of the device under various operating conditions. Figure 5-5 represents such characteristics with the parameters specified in Figure 5-4 under various voltage differences between the emitter and the collector (𝑉𝑎𝑝). 𝜓𝑀 represents the maximum motive and 𝐽𝑡𝑜𝑡 is the total charge flux from the emitter to the collector. As it can be observed from this figure, the overall effect of the introduction of ions is the reduction of space charge. The total current in the presence of ions is about 10 times larger than with no ions in small biases. 103 Figure 5-5. The current density-voltage (𝑉𝑎𝑝) characteristics of the TEC with the parameters specified in Figure 5-2 under various voltage differences between the emitter and the collector. 𝛽 = 0 corresponds to the case where no ions are present. 𝜓𝑀 represents the maximum motive and 𝐽𝑡𝑜𝑡 is the total current density from the emitter to the collector. Summary In summary, the physics of the operation of TECs in low pressure-plasmas was developed. It was observed that, depending on the ratio of the flux density of the ions to the electrons, a wide array of motive profiles could arise. When this ratio is small, the effect of the ions is the total reduction of the space charge. However, as this ratio is enhanced, plasma oscillations could occur. We also demonstrated the extendibility of our model by including electron and ion emissions from both the emitter and the collector. 104 Chapter 6- Contributions and future work The works presented in this thesis provide us with the necessary tools to analyze the behaviour of thermionic converters under several important conditions, including the presence of a grid and or a low-pressure plasma in the device. As explained in the introduction chapter, the pre-existing models could not capture many of the significant demands of nascent, nanomaterial-based TECs. Some of these requirements stem from the presence of a high temperature gradient on the cathode, the need to solve the reverse problem to calculate the device parameters from experimental current-voltage measurements, splicing errors, or operation in the presence of a grid or a low-pressure plasma. The models developed here address all these situations. Contributions In Chapter 2, the pre-existing model of the analysis of TECs that is based on finding an analytical solution to the Vlasov-Poisson system was substantially improved by: o Solving the integrals for a much wider range of parameters o Proposing a self-consistent method to calculate the output characteristics rather than using the approximations that are suggested in [14]. This improvement allowed us to solve the reverse problem and calculate the internal parameters, e.g. workfunctions, inter-electrode distances, effective emission area and temperatures of the thermionic converters. It was further demonstrated that nascent thermionic devices, though operating on the same overall principles, could have some distinct features such as spatially varying workfunctions and the presence of temperature gradients that become pronounced in the presence of localized heating, which require more elaborate simulation methods. Due to the high stability of the proposed method, a wide range of parameters can be solved for and, therefore, the method can address the numerical difficulties associated with emergent thermionic devices. As an example of its application, this method was used to extract the parameters of a light induced thermionic conversion device based on a MWCNT forest (Appendix A). In Chapter 3, we proposed an iterative and self-consistent approach with the following advantages: o Completely avoids the splicing errors that arise in the analytical approach that uses different sets of equations at various operation regimes. o By introducing asymptotic expansion to deal with the ill-behaving functions, the range of applicability of our approach is significantly increased in the presence of an analytical solution. 105 o This method is applicable in cases where a closed-form solution to the Vlasov equation cannot be found. o This model is applicable in cases that the analytical solution can easily become intractable, such as the presence of back-emission from the collector. o A particle tracing approach can be employed to replace the Vlasov equation. This developed model was used as the launching pad for the following Chapters to attack more elaborate problems such as the presence of a grid or the operation of the TEC in the presence of low-pressure plasmas. In Chapter 4, we developed the physics of thermionic converters in the presence of a grid in their inter-electrode region. We demonstrated that: o A momentum gap could arise in the phase-space of the electrons in cases where a potential well is formed in the inter-electrode region. The emergence of this momentum gap was attributed to the reflection of electrons from the edges of the well. o There is a subtle interplay between the concentration of the electrons and their velocity in the presence of these momentum gaps. In other words, the performance of a TEC cannot be captured by the same equations that are applied to a TEC in the absence of grids. This model was employed to calculate the output characteristics of TECs under varying operation conditions, and it was shown that the current output density of a TEC can be improved by 3 orders of magnitude by a proper design of the grid. In Chapter 5, the physics of thermionic converters in the presence of low pressure plasmas was developed. It was shown that: o Depending on the concentrations of the electrons and ions, momentum gaps could arise in the phase space of both electrons and ions. o Depending on the ratio of the saturation current density of the ions to electrons, 𝛽, different scenarios could arise: At low values of 𝛽, such that momentum gaps do not arise in the phase-space of either electrons or ions, the motive is simply lowered while maintaining the same overall profile as in the case that no ions are present. If the value of 𝛽 is high enough such that momentum gaps arise in the phase-space of ions and not that of electrons, the motive will be negative at some points, while positive at others, and the overall effect is the total reduction of the maximum motive in the inter-electrode region. At higher values of 𝛽, momentum gaps could arise in the phase space of both electrons and ions. In such cases, a common solution to the steady state time-independent Vlasov 106 equation and the Poisson equation is non-existent. Therefore, the system has plasma oscillations. We provided arguments on estimating the amplitude of these oscillations and also their periods. Future work The works of this thesis lay the foundation for further development of emergent TECs under varying operational conditions. Some of the next steps that come to mind are: 1. Solving the Schrodinger equation or using a non-equilibrium Green’s function formalism to subsume the regime where quantum tunneling becomes important at nanoscale inter-electrode distances. To this end, the Vlasov equation is replaced with the Schrodinger equation and the potential energy term is calculated in a mean-field fashion from the Poisson equation. The classical results of this thesis were obtained based on the mean-field approach and the range of their applicability was discussed in Chapter 4. The mean-field approach, in the context of quantum mechanics, neglects two quantum effects: the non-commutativity of the kinetic and potential energy operators that enter the Hamiltonian and the lack of antisymmetry in the wave function. Using Trotter’s expansion of the kinetic and potential energy operators in the density matrix [94] [95] [96], it can be shown that the errors caused by non-commutativity of the kinetic and potential energy operators in calculating thermodynamic potentials (such as free energy) are on the order of 𝜆2𝑘𝐵𝑇< (𝜕𝑉𝜕𝑥)2>𝑐𝑙 , where 𝜆 is the thermal de Broglie wavelength of the electrons. < (𝜕𝑉𝜕𝑥)2>𝑐𝑙 represents the expectation value of the square of the derivative of the potential energy with respect to distance (a classical calculation)—we surmise that this is proportional to 1/𝐿𝐷2 , where 𝐿𝐷 is the Debye length in our problem. Therefore, when 𝜆 ≪ 𝐿𝐷, the non-commutativity of the kinetic energy can be neglected. It can be seen that for 𝜆 =ℎ√2 𝜋 𝑚 𝑘𝐵 𝑇 to be on the order of 𝐿𝐷 = (𝜖02𝑘𝐵32 𝜋𝑚 𝑒2)14 𝑇𝐸3/4𝐽1/2, at a temperature of 2,000 𝐾, workfunctions smaller than ~1 𝑒𝑉—which are extremely low—would be required. The other effect is due to neglecting the antisymmetric nature of the wave function. It is noted that the partition function calculated from the microcanonical ensemble already has an “ad hoc” factor to account for the indistinguishability of the electrons; otherwise, unrealistic effects such as lack of extensivity of entropy (Gibbs Paradox) could arise [50] [97]. The values of corrections can be obtained by calculating the density matrix (and hence the partition function) in states that are made antisymmetrized by application of the exchange operator, 𝑃. It makes sense that 107 these effects become important when the average distance between the electrons is smaller than 𝜆. A rigorous proof can be found in Appendix B of ref. [95]. Interestingly, this corresponds to workfunctions smaller than ~ 0.2 𝑒𝑉. Considering that state-of-the-art workfunctions at the high temperature of 2000 𝐾 are on the order of ~ 2 𝑒𝑉, the methods developed in this thesis can be expanded to include the quantum range without the need to solve the many-body Schrodinger equation in the foreseeable future. It is noted that the Fermi temperature corresponding to electrons in TECs is extremely small (due to low electron concentrations), i.e. ~ 6 𝐾 for a workfunction of 1.5 𝑒𝑉, as opposed to ~ 6 × 105 K in solids. Therefore, it is expected that solids act quantum mechanically (for instance a classical kinetic theory leads to spurious results when calculating the heat capacity of metals) even at very high temperatures, whereas, extremely high densities of electrons are necessary to have appreciable quantum effects (with no classical counterpart) in TECs. 2. Extending the Vlasov equation to cases where a hemi-Maxwellian distribution does not apply. This can be particularly useful when photoemission from the electrodes cannot be neglected. 3. Extending the model to time-dependent cases. This can be used to study the plasma oscillations that were discussed in Chapter 5. 4. Development of a three dimensional particle tracing model. 5. Including magnetic fields, either in particle tracing or in the Vlasov-Poisson system. 6. Experimental validation of the theoretically derived results. Although the models developed in Chapters 2 and 3 have been tested using the experimental results shown in the Appendices, it would be beneficial to corroborate the models of chapters 4 and 5 as well. 108 References [1] A. Khoshaman, H. D. E. Fan, A. Koch, G. . Sawatzky, and A. Nojeh, “Thermionics, Thermoelectrics, and Nanotechnology: New Possibilities for Old Ideas,” IEEE Nanotechnol. Mag., vol. 8, no. 2, pp. 4–15, Jun. 2014. [2] A. H. Khoshaman, A. T. Koch, M. Chang, H. D. E. Fan, M. V. Moghaddam, and A. Nojeh, “Nanostructured Thermionics for Conversion of Light to Electricity: Simultaneous Extraction of Device Parameters,” IEEE Trans. Nanotechnol., vol. 14, no. 4, pp. 624–632, Jul. 2015. [3] A. H. Khoshaman and A. Nojeh, “A self-consistent approach to the analysis of thermionic devices,” J. Appl. Phys., vol. 119, no. 4, p. 44902, Jan. 2016. [4] A. H. Khoshaman and A. Nojeh, “Classical momentum gap for electron transport in vacuum and consequences for space charge in thermionic converters with a grid electrode,” J. Vac. Sci. Technol. B, vol. 34, no. 4, p. 40610, Jul. 2016. [5] Khoshaman, A.H. and Nojeh, A., “Low-pressure plasma-enhanced behavior of thermionic convertors,” Journal of applied physics, submitted. [6] “Publication: Key World Energy Statistics 2015.” [Online]. Available: https://www.iea.org/publications/freepublications/publication/key-world-energy-statistics-2015.html. [Accessed: 14-Mar-2016]. [7] “Laser-induced thermoelectric energy generation using carbon nanotube forests - UBC Library Open Collections.” [Online]. Available: https://open.library.ubc.ca/cIRcle/collections/ubctheses/24/items/1.0167189. [Accessed: 23-Mar-2016]. [8] “Global Exergy Resource Chart - GCEP.” [Online]. Available: http://gcep.stanford.edu/research/exergy/resourcechart.html. [Accessed: 14-Mar-2016]. [9] V. Smil, Power Density: A Key to Understanding Energy Sources and Uses. Cambridge, Massachusetts: The MIT Press, 2015. [10] J. H. Lee, I. Bargatin, T. O. Gwinn, M. Vincent, K. A. Littau, R. Maboudian, Z.-X. Shen, N. A. Melosh, and R. T. Howe, “Microfabricated silicon carbide thermionic energy converter for solar electricity generation,” in 2012 IEEE 25th International Conference on Micro Electro Mechanical Systems (MEMS), 2012, pp. 1261–1264. [11] P. Yaghoobi, M. V. Moghaddam, and A. Nojeh, “‘Heat trap’: Light-induced localized heating and thermionic electron emission from carbon nanotube arrays,” Solid State Commun., vol. 151, no. 17, pp. 1105–1108, Sep. 2011. [12] “A reformulation of thermionic theory for vacuum diodes.” [Online]. Available: http://www.sciencedirect.com/science/article/pii/S0039602802020630. [Accessed: 16-Mar-2016]. [13] N. S. Rasor, “Thermionic energy conversion plasmas,” Plasma Sci. IEEE Trans. On, vol. 19, no. 6, pp. 1191–1208, Dec. 1991. [14] G. N. Hatsopoulos and E. P. Gyftopoulos, Thermionic Energy Conversion, Vol. 1. The MIT Press, 1973. [15] S. Meir, C. Stephanos, T. H. Geballe, and J. Mannhart, “Highly-efficient thermoelectronic conversion of solar energy and heat into electricpower,” J. Renew. Sustain. Energy, vol. 5, no. 4, p. 43127, Jul. 2013. [16] R. M. Williams, Electrical performance of a rhenium-niobium cylindrical thermionic converter /. Washington, D. C. :, 1968. [17] “Power Density Primer: Understanding the Spatial Dimension of the Unfolding Transition to Renewable Electricity Generation (Part IV – New Renewables Electricity Generation) - Master Resource.” [Online]. Available: https://www.masterresource.org/smil-vaclav/smil-density-new-renewables-iv/. [Accessed: 22-Apr-2016]. 109 [18] “Samsung S7.0-171 - 7.000,0 kW - Turbine.” [Online]. Available: http://en.wind-turbine-models.com/turbines/285-samsung-s7.0-171. [Accessed: 01-Aug-2016]. [19] T. Jarvis, “Thermionic conversion-state of the art,” IRE Trans. Mil. Electron., vol. MIL-6, no. 1, pp. 41–45, Jan. 1962. [20] C. B. Vining, “An inconvenient truth about thermoelectrics,” Nat. Mater., vol. 8, no. 2, pp. 83–85, Feb. 2009. [21] “Steam Turbine Efficiency,” Turbinesinfo - All About Turbines, 07-Aug-2011. . [22] W. Xiang and Y. Chen, “Performance improvement of combined cycle power plant based on the optimization of the bottom cycle and heat recuperation,” J. Therm. Sci., vol. 16, no. 1, pp. 84–89. [23] T. M. Tritt and M. A. Subramanian, “Thermoelectric Materials, Phenomena, and Applications: A Bird’s Eye View,” MRS Bull., vol. 31, no. 3, pp. 188–198, Mar. 2006. [24] I. N. Levine, Physical chemistry, 5th ed. Boston: McGraw-Hill, 2002. [25] J. H. Ingold, “Calculation of the Maximum Efficiency of the Thermionic Converter,” J. Appl. Phys., vol. 32, no. 5, pp. 769–772, May 1961. [26] N. S. Rasor, “Thermionic energy conversion plasmas,” IEEE Trans. Plasma Sci., vol. 19, no. 6, pp. 1191–1208, Dec. 1991. [27] J. M. Houston, “Theoretical Efficiency of the Thermionic Energy Converter,” J. Appl. Phys., vol. 30, no. 4, pp. 481–487, Apr. 1959. [28] B. Zhao, H. Hu, A. Yu, D. Perea, and R. C. Haddon, “Synthesis and Characterization of Water Soluble Single-Walled Carbon Nanotube Graft Copolymers,” J. Am. Chem. Soc., vol. 127, no. 22, pp. 8197–8203, Jun. 2005. [29] A. H Khoshaman, H. D. Fan, and A. Nojeh, “Improving the Efficiency of Light Induced Thermionic Energy Conversion of Carbon Nanotubes by Potassium Intercalation,” presented at the 97th Canadian Chemistry Conference and Exhibition, Vancouver, BC. [30] J. R. Smith, G. L. Bilbro, and R. J. Nemanich, “Theory of space charge limited regime of thermionic energy converter with negative electron affinity emitter,” J. Vac. Sci. Technol. B, vol. 27, no. 3, pp. 1132–1141, May 2009. [31] V. S. Robinson, T. S. Fisher, J. A. Michel, and C. M. Lukehart, “Work function reduction of graphitic nanofibers by potassium intercalation,” Appl. Phys. Lett., vol. 87, no. 6, p. 61501, 2005. [32] J. W. Schwede, I. Bargatin, D. C. Riley, B. E. Hardin, S. J. Rosenthal, Y. Sun, F. Schmitt, P. Pianetta, R. T. Howe, Z.-X. Shen, and N. A. Melosh, “Photon-enhanced thermionic emission for solar concentrator systems,” Nat. Mater., vol. 9, no. 9, pp. 762–767, 2010. [33] A. N. Tavkhelidze, “Nanostructured electrodes for thermionic and thermotunnel devices,” J. Appl. Phys., vol. 108, no. 4, p. 44313, 2010. [34] A. Shakouri, Z. Bian, R. Singh, Y. Zhang, D. Vashaee, T. E. Humphrey, H. Schmidt, J. M. Zide, G. Zeng, and J.-H. Bahk, “Solid-State and Vacuum Thermionic Energy Conversion,” CALIFORNIA UNIV SANTA CRUZ SCHOOL OF ENGINEERING, CALIFORNIA UNIV SANTA CRUZ SCHOOL OF ENGINEERING, 2005. [35] D. B. Go, “Theoretical analysis of ion-enhanced thermionic emission for low-temperature, non-equilibrium gas discharges,” J. Phys. Appl. Phys., vol. 46, no. 3, p. 35202, Jan. 2013. [36] B. Y. Moyzhes and T. H. Geballe, “The thermionic energy converter as a topping cycle for more efficient heat engines—new triode designs with a longitudinal magnetic field,” J. Phys. Appl. Phys., vol. 38, no. 5, pp. 782–786, Mar. 2005. [37] S. Basu, Z. M. Zhang, and C. J. Fu, “Review of near-field thermal radiation and its application to energy conversion,” Int. J. Energy Res., vol. 33, no. 13, pp. 1203–1232, 2009. [38] J.-H. Lee, I. Bargatin, N. A. Melosh, and R. T. Howe, “Optimal emitter-collector gap for thermionic energy converters,” Appl. Phys. Lett., vol. 100, no. 17, p. 173904, 2012. [39] J. I. Lee, Y. H. Jeong, H.-C. No, R. Hannebauer, and S.-K. Yoo, “Size effect of nanometer vacuum gap thermionic power conversion device with CsI coated graphite electrodes,” Appl. Phys. Lett., vol. 95, no. 22, p. 223107, 2009. 110 [40] R. Yang, A. Narayanaswamy, and G. Chen, “Surface-Plasmon Coupled Nonequilibrium Thermoelectric Refrigerators and Power Generators,” J. Comput. Theor. Nanosci., vol. 2, no. 1, pp. 75–87, Mar. 2005. [41] Z. Taliashvili, A. Tavkhelidze, L. Jangidze, and Y. Blagidze, “Vacuum nanogap formation in multilayer structures by an adhesion-controlled process,” Thin Solid Films, vol. 542, pp. 399–403, Sep. 2013. [42] T. Zeng, “Thermionic-tunneling multilayer nanostructures for power generation,” Appl. Phys. Lett., vol. 88, no. 15, p. 153104, 2006. [43] S. F. Adams, “Solar Thermionic Space Power Technology Testing: A Historical Perspective,” AIP Conf. Proc., vol. 813, no. 1, pp. 590–597, Jan. 2006. [44] K. L. Jensen, “General formulation of thermal, field, and photoinduced electron emission,” J. Appl. Phys., vol. 102, no. 2, p. 24911, 2007. [45] R. H. Fowler, “The Analysis of Photoelectric Sensitivity Curves for Clean Metals at Various Temperatures,” Phys. Rev., vol. 38, no. 1, pp. 45–56, Jul. 1931. [46] L. A. DuBridge, “Theory of the Energy Distribution of Photoelectrons,” Phys. Rev., vol. 43, no. 9, pp. 727–741, May 1933. [47] K. Nakayama, K. Tanabe, and H. A. Atwater, “Plasmonic nanoparticle enhanced light absorption in GaAs solar cells,” Appl. Phys. Lett., vol. 93, no. 12, p. 121904, 2008. [48] P. Yaghoobi, M. Vahdani Moghaddam, and A. Nojeh, “Solar electron source and thermionic solar cell,” AIP Adv., vol. 2, no. 4, pp. 42139-42139–12, Nov. 2012. [49] P. Yaghoobi, “Laser-induced electron emission from arrays of carbon nanotubes,” University of British Columbia, 2012. [50] M. Tuckerman, Statistical Mechanics: Theory and Molecular Simulation. Oxford ; New York: Oxford University Press, 2010. [51] R. C. Tolman, The principles of statistical mechanics. London: Oxford University Press, 1938. [52] F. F. Chen, Introduction to plasma physics and controlled fusion. Volume 1, Plasma physics, 2nd edition. New York: Springer, 2006. [53] D. L. Pulfrey, Understanding Modern Transistors and Diodes, 1 edition. Cambridge ; New York: Cambridge University Press, 2010. [54] J. A. Bittencourt, Fundamentals of Plasma Physics. New York, NY: Springer New York, 2004. [55] H. Goldstein, C. P. Poole, and J. L. Safko, Classical mechanics, 3rd ed. San Francisco: Addison Wesley, 2002. [56] I. Langmuir, “The Effect of Space Charge and Initial Velocities on the Potential Distribution and Thermionic Current between Parallel Plane Electrodes,” Phys. Rev., vol. 21, no. 4, pp. 419–435, Apr. 1923. [57] G. N. Hatsopoulos and E. P. Gyftopoulos, Thermionic Energy Conversion - Vol. 2: Theory, Technology, and Application. The MIT Press, 1979. [58] D. J. Griffiths, Introduction to Electrodynamics, 4 edition. Boston: Pearson, 2012. [59] A. Zangwill, Modern Electrodynamics, 1 edition. Cambridge: Cambridge University Press, 2012. [60] J. D. Jackson, Classical Electrodynamics Third Edition, 3 edition. New York: Wiley, 1998. [61] J. R. Smith, G. L. Bilbro, and R. J. Nemanich, “Considerations for a high-performance thermionic energy conversion device based on a negative electron affinity emitter,” Phys. Rev. B, vol. 76, no. 24, p. 245327, Dec. 2007. [62] G. L. B. Joshua Ryan Smith, “Theory of space charge limited regime of thermionic energy converter with negative electron affinity emitter,” J. Vac. Sci. Amp Technol. B Microelectron. Nanometer Struct., vol. 27, no. 3, 2009. [63] J. R. Smith, “Increasing the Efficiency of a Thermionic Engine Using a Negative Electron Affinity Collector,” J. Appl. Phys., vol. 114, no. 16, p. 164514, 2013. [64] C. Stephanos, “Thermoelectronic Power Generation from Solar Radiation and Heat,” Universität Augsburg, 2012. 111 [65] F. S. Lo, P. S. Lu, B. Ragan-Kelley, A. Minnich, T. H. Lee, M. C. Lin, and J. P. Verboncoeur, “Modeling a thermionic energy converter using finite-difference time-domain particle-in-cell simulations,” Phys. Plasmas 1994-Present, vol. 21, no. 2, p. 23510, Feb. 2014. [66] S. Srisonphan, M. Kim, and H. K. Kim, “Space charge neutralization by electron-transparent suspended graphene,” Sci. Rep., vol. 4, Jan. 2014. [67] S. Suzuki, Y. Watanabe, Y. Homma, S. Fukuba, S. Heun, and A. Locatelli, “Work functions of individual single-walled carbon nanotubes,” Appl. Phys. Lett., vol. 85, no. 1, pp. 127–129, Jul. 2004. [68] M. Grujicic, G. Cao, and B. Gersten, “Enhancement of field emission in carbon nanotubes through adsorption of polar molecules,” Appl. Surf. Sci., vol. 206, no. 1–4, pp. 167–177, Feb. 2003. [69] A. Maiti, J. Andzelm, N. Tanpipat, and P. von Allmen, “Effect of Adsorbates on Field Emission from Carbon Nanotubes,” Phys. Rev. Lett., vol. 87, no. 15, p. 155502, 2001. [70] G. Zhou, W. Duan, and B. Gu, “Electronic Structure and Field-Emission Characteristics of Open-Ended Single-Walled Carbon Nanotubes,” Phys. Rev. Lett., vol. 87, no. 9, p. 95504, Aug. 2001. [71] D. Deutsch, The beginning of infinity: explanations that transform the world. New York: Viking, 2011. [72] D. Deutsch, The fabric of reality, 1st. ed. Harmondsworth: Allan Lane, 1997. [73] A. F. Dugan, “Contribution of Anode Emission to Space Charge in Thermionic Power Converters,” J. Appl. Phys., vol. 31, no. 8, pp. 1397–1400, Aug. 1960. [74] J. A. Michel, V. S. Robinson, L. Yang, S. Sambandam, W. Lu, T. Westover, T. S. Fisher, and C. M. Lukehart, “Synthesis and characterization of potassium metal/graphitic carbon nanofiber intercalates,” J. Nanosci. Nanotechnol., vol. 8, no. 4, pp. 1942–1950, Apr. 2008. [75] T. L. Westover, A. D. Franklin, B. A. Cola, T. S. Fisher, and R. G. Reifenberger, “Photo- and thermionic emission from potassium-intercalated carbon nanotube arrays,” J. Vac. Sci. Technol. B Microelectron. Nanometer Struct., vol. 28, no. 2, pp. 423–434, 2010. [76] A. H. Khoshaman, H. D. E. Fan, A. T. Koch, N. H. Leung, and A. Nojeh, “Localized light induced thermionic emission from intercalated carbon nanotube forests,” in Vacuum Nanoelectronics Conference (IVNC), 2014 27th International, 2014, pp. 59–60. [77] W. M. Haynes, Ed., CRC Handbook of Chemistry and Physics, 93rd Edition, 93rd ed. CRC Press, 2012. [78] P. L. Auer, “Potential Distributions in a Low‐Pressure Thermionic Converter,” J. Appl. Phys., vol. 31, no. 12, pp. 2096–2103, Dec. 1960. [79] Y. Li, R. Tirumala, P. Rumbach, and D. B. Go, “The coupling of ion-enhanced field emission and the discharge during microscale breakdown at moderately high pressures,” IEEE Trans. Plasma Sci., vol. 41, no. 1, pp. 24–35, 2013. [80] D. B. Go, T. S. Fisher, S. V. Garimella, and V. Bahadur, “Planar microscale ionization devices in atmospheric air with diamond-based electrodes,” Plasma Sources Sci. Technol., vol. 18, no. 3, p. 35004, 2009. [81] R. Tirumala and D. B. Go, “Multi-electrode assisted corona discharge for electrohydrodynamic flow generation in narrow channels,” IEEE Trans. Dielectr. Electr. Insul., vol. 18, no. 6, pp. 1854–1863, 2011. [82] C. R. Crowell, “The Richardson constant for thermionic emission in Schottky barrier diodes,” Solid-State Electron., vol. 8, no. 4, pp. 395–399, Apr. 1965. [83] “Principles Plasma Physics Engineers And Scientists | Electromagnetics | Cambridge University Press.” [Online]. Available: http://www.cambridge.org/us/academic/subjects/engineering/electromagnetics/principles-plasma-physics-engineers-and-scientists. [Accessed: 28-Jul-2014]. [84] M. E. Kiziroglou, X. Li, A. A. Zhukov, P. A. J. de Groot, and C. H. de Groot, “Thermionic field emission at electrodeposited Ni–Si Schottky barriers,” Solid-State Electron., vol. 52, no. 7, pp. 1032–1038, Jul. 2008. [85] J. Orloff, Ed., Handbook of charged particle optics, 2nd ed. Boca Raton: CRC Press/Taylor & Francis, 2009. 112 [86] “M. Chang, M. Vahdani Moghaddam, A. Khoshaman, M.S. Mohammed Ali, M. Dahmardeh, K. Takahata, A. Nojeh, ‘High Temperature Gradient in a Conductor: Carbon Nanotube Forest under the ’Heat Trap‘ Condition,’ 57th Int’l Conf. Electron, Ion, Photon Beam Technology and Nanofabrication (EIPBN 2013), Nashville, USA, 2013, http://eipbn.omnibooksonline.com/data/papers/2013/P11-08.pdf#page=1 (2pp).” [87] S.-J. Liang and L. K. Ang, “Electron Thermionic Emission from Graphene and a Thermionic Energy Converter,” Phys. Rev. Appl., vol. 3, no. 1, p. 14002, Jan. 2015. [88] G. B. Arfken and H. J. Weber, Mathematical Methods for Physicists, Seventh Edition: A Comprehensive Guide, 7 edition. Amsterdam ; Boston: Academic Press, 2012. [89] R. D. Skeel, “Variable step size destabilizes the Störmer/leapfrog/verlet method,” BIT Numer. Math., vol. 33, no. 1, pp. 172–175, Mar. 1993. [90] C.K. Birdsall, A.B. Langdon, Plasma Physics via Computer Simulation, 1 edition. CRC Press, 2004. [91] L. D. Landau and L. D. Landau, Mechanics. Oxford: Pergamon, 1960. [92] G. L. B. Joshua Ryan Smith, “Theory of space charge limited regime of thermionic energy converter with negative electron affinity emitter,” J. Vac. Sci. Amp Technol. B Microelectron. Nanometer Struct., vol. 27, 2009. [93] W. M. McKeeman, “Algorithm 145: Adaptive Numerical Integration by Simpson’s Rule,” Commun ACM, vol. 5, no. 12, p. 604–, Dec. 1962. [94] M. Suzuki, Quantum Monte Carlo Methods in Equilibrium and Nonequilibrium Systems Proceedings of the Ninth Taniguchi International Symposium, Susono, Japan, November 14-18, 1986. S.l.: Springer Berlin Heidelberg, 1987. [95] F. Schwabl, Statistical mechanics. Berlin ; New York: Springer, 2002. [96] F. Schwabl, Quantum mechanics, 3rd ed. Berlin ; New York: Springer, 2002. [97] C. Kittel, Elementary statistical physics. --. New York: Wiley, 1958. [98] P. Liu, Y. Wei, K. Jiang, Q. Sun, X. Zhang, S. Fan, S. Zhang, C. Ning, and J. Deng, “Thermionic emission and work function of multiwalled carbon nanotube yarns,” Phys. Rev. B, vol. 73, no. 23, Jun. 2006. [99] P. Liu, Q. Sun, F. Zhu, K. Liu, K. Jiang, L. Liu, Q. Li, and S. Fan, “Measuring the Work Function of Carbon Nanotubes with Thermionic Method,” Nano Lett., vol. 8, no. 2, pp. 647–651, Feb. 2008. 113 Appendix A: Light-induced thermionic converter experimental results and comparisons with theory6 As an experimental case for the model developed in Chapter 2, a carbon nanotube forest is used as the emitter in a light based thermionic energy converter and locally heated to thermionic emission temperatures using a 50-mW-focused laser beam. The current-voltage characteristics were measured and used to solve the reverse problem to obtain the internal parameters of the device, which were shown to be consistent with the values obtained using other methods. A.1. Experimental set-up Millimeter-long MWCNT forests were grown in a home-built chemical vapor deposition reaction chamber over a highly p-doped silicon wafer coated with successive layers of 10 nm of Alumina and 1 nm of Iron. Some of the catalysts were patterned to form cubically shaped forests with a base width of 1 mm separated by a distance of ~500 μm. This pattern facilitated the characterization of the illuminated spot on the sidewalls of the CNT forests. CNTs were grown on the catalyst by following a recipe similar to the one used in [48]. The following experimental set-up (Figure A-1(a)) was adopted for the experiments. The chips containing CNT forests were mounted on a sample holder and placed in a high-vacuum chamber (10-5 – 10-6 torr). The vacuum chamber was equipped with Pirani and cold cathode gauges that allowed the pressure to be monitored within a wide range (103 – 10-3 torr and 10-2 – 10-10 torr respectively). A precision leak valve was used to control the pressure inside the chamber. The CNT sample was used as the cathode and was illuminated on its sidewall through the anode. Indium thin oxide (ITO) glass slide with 84 % transparency or an 85% transparent stainless steel mesh were used as the anode. The gap width was set to less than 1 mm. The cathode was surrounded by an open-ended aluminum cube, kept at the same potential as the anode. 6 The results of this Appendix are reproduced with minor modifications from [2]. 114 The laser was focused on the sidewalls of the CNT forest by means of a lens with a 15-cm focal length. Two manual high precision microactuators were used to move laterally along the width of the forest, allowing one to determine the location of the illuminated spot for characterization purposes. The lens was placed approximately 15 cm away from the anode through a viewport and its exact position was controlled by a high precision motorized microactuator (with a precision of 200 nm) controlled by a computer. Several different types of lasers, including three lasers from Laserglow Technologies (Electra Pro-40 (405 nm), Scorpius-500 (1064 nm) and Hercules (532 nm)) and a 300- mW, 532-nm Green DPSS Laser from Alltech solutions were used. The set-up included a broadband non-polarizing hybrid cube beamsplitter allowing the alignment of two laser beams to impinge on the sidewall of the forest at the same spot. Two optical attenuators and a power reader were used to adjust the power to the desired value. A charged coupled device (CCD) camera was used to measure the area of the laser spot on the forest (Figure A-1(b)). The I-V characteristics were measured by sweeping from negative (retarding) to positive (collection) voltages using Keithley 6430 and 6517a electrometers. Since the voltage is swept over a wide range in our experiments, the magnitude of current varies dramatically each time (typically from less than 1 pA to ~10 μA for a voltage sweep in the range of -4 – 10V at 𝑇𝐸 = 2,300 𝐾). At positive biases higher than the saturation voltage of the LITEC device (typically 𝑉 = 30 𝑉 to ensure that the saturation regime is properly covered), the current is on the order of tens of 𝜇𝐴s. 115 Figure A-1. (a) The schematic of the experimental set-up for I-V measurements of a CNT-based LITEC. (b) The incandescent spot on the side-wall of the CNT forest as imaged by a CCD camera. (c) The internal parameters of the LIETC device are obtained by solving the reverse problem and fitting to the experimental results as depicted in the legend of the plot. A.2. Extracting the experimental parameters As explained in the theoretical section, each regime of operation follows a certain set of equations. Each internal parameter of the LITEC device affects the I-V characteristics in a unique way. For instance, increasing the gap distance (from Figure 2-4) increases the width of the space charge region while maintaining the same level of saturation current. Raising the emitter’s temperature shifts the entire I-V curve upwards (increasing current at an exponential rate), while maintaining the intersection of the two tangent lines to the characteristics in the retarding and the saturation modes at a fixed voltage point, 𝑉𝑚 =(𝜙𝐶 − 𝜙𝐸)/𝑒. Changing 𝜙𝐸 affects the saturation point voltage and current while altering the thermionic emission area’s radius (𝑟), shifting the entire curve upward or downward. To distinguish between the 116 alterations caused by 𝑇𝐸 or 𝑟, it is noted that the former has a more palpable effect in negative biases, while the latter shifts the curve uniformly. Although it seems that each set of parameters will lead to only a unique I-V profile, solving for each of the parameters is a challenging problem. In order to solve this reverse problem, a graphical user interface (GUI) was developed to calculate the I-V characteristics based on the proposed procedure, allowing one to select different temperature profiles. Each of the internal parameters are adjusted while their impact on the I-V profile is studied in real time. In order to test the uniqueness, a set of reasonable values of internal parameters were chosen in random, the I-V curves were generated based on the proposed procedure in the previous section, and the results were used to solve the reverse problem. When the temperature was assumed to be uniform, the reverse problem resulted in unique estimations of the internal parameters, with ±0.1 eV of accuracy for workfunctions, ~± 50 K for 𝑇𝐸 and ~± 20 μm for the spot-size radius. The estimation of the inter-electrode distance depends on how deep the device is within the space charge region. At relatively low d values, where 𝑉𝐶 is close to 𝑉𝑠𝑎𝑡, the length of the space charge region is negligible and, therefore, the impact of the gap width on the I-V profile will not be noticeable. Consequently, only an upper boundary for d can be determined in such cases. Figure A-1(c) shows the experimental I-V characteristics from illumination of the sample (with a stainless steel mesh as the anode) by 40 mW of a 405-nm laser fitted to simulation curves identified by a fixed-gradient temperature profile. The effective hot spot radius (121 μm) is comparable with the image taken by a CCD camera (~150 μm). The estimated temperature (2,150 K) when assuming a linear temperature distribution along the hot spot matches closely with the results that our group previously obtained from fitting the incandescence spectra to black body radiation [11], and the estimated workfunction of CNTs (4.7 eV) is consistent with the values obtained by other methods [98]. It is noted that the sidewall of the CNT forest comprises of ~ 105 individual CNTs and since the minimum spatial range we are interested in is on the order of ~ 10 𝜇𝑚, it already contains many individual CNTs. Therefore, the granularity of the heated area can be neglected and a continuum analysis would suffice for our purposes. Since the modelling presented in this work is based on a one dimensional potential profile with infinitely large electrodes, it requires some justification to use this model to extract the parameters of a CNT-based LITEC. Here, we show that our 1-D model can effectively capture the experimental conditions by a relatively small error. These steps are followed to calculate the error due to a 1-D model: 117 1. The electric field in the inter-electrode distance in front of every ring is due to two sources: the external applied bias, and the field due to the presence of electrons, or the space charge field. Therefore, in the absence of the space charge effects, the field in front of every ring is equal. 2. As mentioned, the space charge is identified by the critical point voltage, 𝑉𝐶 and the saturation point voltage, 𝑉𝑠𝑎𝑡. Under the foregoing experimental conditions, only a small number of the rings inhabit the space charge dominated regime for a mere fraction of their entire IV characteristics. For instance, in the case of the lowest-temperature ring (1500 𝐾), the 𝑉𝐶 and 𝑉𝑠𝑎𝑡 coincide at 1.2001 𝑉 (as calculated using the proposed method in the Chapter), signifying that this ring does not enter the space charge regime at all. Similar calculations reveal that the device enters the space charge effect regime of operation at temperatures higher than 1900 𝐾, with 𝑉𝐶 = 1.0654 𝑉 and 𝑉𝑠𝑎𝑡 =1.3651 𝑉. Therefore, the 1-D approximation is fully valid for about 3/4 of the rings. For the remaining 1/4 of the rings, the approximation is not valid for a very short portion of applied biases (for instance 0.3 𝑉 on 40 𝑉 range of voltage sweep). 3. The width of the space charge regime is narrower in a CNT-based LITE device as opposed to the traditional case of a hot electrode in an infinitely long parallel plate capacitor setting. The reason is that the charge density in the transverse direction reduces with distance from the axis, since electrons only are fired from the finite size hot spot. Therefore, for instance, we would expect that the length of the space charge limited regime at 1900 𝐾 to be less than 0.3 𝑉, an estimation of how much lower it should be is presented in the next point. 4. The estimation of the maximum change can be obtained as follows: in our experimental case, a hot spot with area 𝐴(0) = 𝜋𝑅2 and average electron density 𝑛(0), expands to an approximate area of 𝐴(𝑥) = 𝜋(𝑅 + 𝑣𝑡ℎ × 𝑡𝑡𝑜𝑓)2 with its electron density diluted to 𝑛2(𝑥) = 𝑛1(𝑥) × 𝐴(0)/𝐴(𝑥), where 𝑣𝑡ℎ = (2 𝑘 𝑇/𝑚)1/2, 𝑡 is the average time of flight of the electrons (𝑡 =𝑥<𝑣𝑒𝑥> ), < 𝑣𝑒𝑥 > is the average velocity of electrons along the x direction and 𝑛1(𝑥) is the density of electrons at position x in the ideal case that can be calculated using equation 2 and Figure 3. On the other hand, in the ideal case, a portion of the hot spot with area 𝐴(0) maintains its surface area due to the charge migration from the adjacent circles and its electron density, 𝑛1(𝑥) varies along the gap distance according to equation 2. The potential of a uniform circular surface charge distribution, 𝜎(𝑥) =𝑛(𝑥) 𝑑𝑥 at a distance x from the center of the circle can be easily calculated as 𝑉 =𝜎(𝑥)2𝜖0(√𝑅2 + 𝑥2 − 𝑥). Therefore, the fractional change in the potential barrier at the emitter due to dilution of electrons at x outside the collector can be calculated as 𝑛2(𝑥) 𝑑𝑥 2𝜖0(√(𝑅 + 𝑣𝑡ℎ × 𝑡)2 + 𝑥2 − 𝑥). This potential was numerically integrated along the gap 118 distance to calculate the overall change the emitter’s potential due to dilution of the electron density due to finiteness of the heated area. The same procedure was employed to calculate the generated potential at the emitter in the ideal case. This led to a maximum introduced error of 15 % in the resulting potential occurring at the bias point, 𝑉 = 1.2 𝑉 at 1500 𝐾. Similar calculations including the ring topology with inner and outer radii of 𝑅𝑎 and 𝑅𝑏 respectively, led to a maximum of 12 % of error in the ring’s own potential at the emitter and 7 % of error in the potential of the adjacent ring. 5. We can infer that the width of space charge-charge limited regime of operation is roughly 12 +7 = 19 % of the width predicted by the ideal model (for instance instead of 0.3 𝑉, it is 0.24 𝑉 wide) . Based on these observations, the inconsistency between the experimental and simulated points occurs in only 2 % of the data points. Simultaneous measurement of all of these parameters that pertain to a certain experiment in LITEC has not been possible in the past. For instance, in order to measure the workfunction, a common method such as ultraviolet spectroscopy (UPS) can be employed. However, considering that the exact workfunction depends on the local morphology of the carbon nanotubes and the temperature [99], it would be extremely challenging, if not impossible, to have the illuminated spot on the CNTs with exactly the same amount of intensity as that used in the thermionic conversion experiment. Besides the difficulty in determining the actual workfunction, the temperature may vary widely depending on the local density of the carbon nanotube forest. For instance, the Richardson-Dushman constant calculated from the thermionic emission data from [68] were in discrepancy with values obtained from metals. This issue was attributed to the extra surface roughness and different surface area from CNTs compared to metals. Furthermore, measurement of the hot-spot area depends on the resolution as well as the bandwidth of the CCD camera used for the measurements. Therefore, the area measured by a CCD camera might be different from the effective area contributing to thermionic emission. However, using I-V curves and following the proposed strategy provides a consistent method for determination of the effective area of the hot spot. Moreover, the effective gap width between the electrodes, though being on the same order as its equivalent physical value, can be different due to the collection efficiency of the electrons. The collection efficiency depends on the configuration of the electrodes and the geometry of the structure. As a result, the value of the gap width obtained from fitting to simulation curves can be interpreted as the effective gap width. The magnitude of this parameter is proportional to the physical inter-electrode distance divided by the collection coefficient. 119 Appendix B: Light induced thermionic energy conversion using in-situ intercalated carbon nanotube forests7 Here we report on the influence of potassium intercalation on the workfunction of CNTs under LITE conditions. Due to the extremely reactive nature of potassium, all experiments were performed under an inert argon environment inside a glove-box. For the ex-situ experiments, potassium was deposited onto the sidewalls of nanotube forests in a custom-built vacuum chamber while being heated to 423 K for a period of 48 hours and pumped to a pressure of 10−3𝑇𝑜𝑟𝑟 The remaining potassium metal retained its metallic shininess, indicating the successful removal of oxygen gas from the reaction chamber. Figure B-1. The sample holder apparatus used for the in-situ intercalation of carbon nanotubes under LITE conditions. Potassium is placed inside a custom-built electrically insulating crucible which is mounted on a resistive heater. The potassium is evaporated and deposited on the surface and side-walls of multiwalled 7This chapter is reproduced with minor modifications from [76] and [29]. 120 carbon nanotubes. The entire set-up is placed inside a custom-built vacuum chamber. The output current-voltage characteristics are measured using a Keithley electrometer. The intercalated samples inside the glove-box were then mounted as the cathode inside a high-vacuum chamber, where an ITO electrode served as the anode placed about 100 𝜇𝑚 away from the emitter. The Chamber was pumped down to a pressure of 10−5 𝑇𝑜𝑟𝑟. The intercalated sidewall of the CNT forest was illuminated by various powers of 402 and 532 𝑛𝑚 lasers. The radius of the laser beam on the side-wall was ~ 130 𝜇𝑚. A Keithley electrometer was employed to measure the current-voltage (I-V) characteristics via a sweep from a range of negative (retarding) to positive (collection) voltages (Figure B-1). For the in-situ experiment, potassium metal was placed inside a ceramic crucible on top of a resistive heater inside the glove-box. The CNT forest was mounted directly on top of the potassium reservoir. Opposite to the CNT forest and at an approximate distance of 1 𝑚𝑚, a stainless steel mesh with 85 % transparency was mounted as the anode (Figure B-1).. The potassium reservoir was heated resistively by passing various currents for a duration of 72 hours. The I-V characteristics of the system were measured while the CNT forest was illuminated by various powers of incident light. B.1. Results and discussion Figure B-2 depicts the I-V characteristics of the bare CNT forest sample before it was intercalated. Based on the method that we proposed in [2], the workfunction of the bare CNT forest was calculated as 4.6 eV. 121 Figure B-2. The current-voltage characteristic of the MWCNT sample before the intercalation. The value of the workfunction of the CNTs was calculated based on our previous work [2] to be 4.6 𝑒𝑉. The input laser power for this experiment is 25 𝑚𝑊. Figure B-3 reveals the emission current as a function of time after the intercalation (for an ex-situ condition, i.e., when the resistive heater is turned off and the potassium atoms on the surface are not replenished). The amount of workfunction reduction was calculated from the difference between the initial current and the steady state current based on the Richardson-Dushman equation. The average observed workfunction reduction was estimated to be in the range of 0.5 eV to 1.1 eV, depending on the incident power, and therefore the temperature of the hot spot. A positive bias of 25 V was applied to the ITO electrode to ensure complete extraction of all the electrons that are thermionically emitted from the CNTs. The reduction of current was attributed to de-intercalation of potassium from the samples, since infliction of damage to CNTs at such low currents is unlikely. In the ex-situ case, to ensure that the current reduction is due to potassium de-intercalation, in the beginning high laser input power (100 𝑚𝑊) was used and the current was allowed to stabilize. Subsequently, the laser was turned off and the heater was turned on for an hour to induce evaporation of the potassium atoms on the surface of CNTs. Figure B-3. The time evolution of current between the collector and the emitter in an intercalated sample under ex-situ conditions. The reason behind using an in-situ intercalation of CNTs is to circumvent the de-intercalation problems that can arise in an ex-situ experiment. Having potassium vapour from a heated reservoir leads to absorption/desorption equilibrium in steady state. (Moreover, the space charge effect can be mitigated in the presence of positive ions, which is another significant advantage.) This condition was achieved in the in-situ intercalation experiment and the corresponding I-V characteristics are shown in Figure B-4. Based 122 on the data shown on this figure, we can calculate the new value of workfunction and see that it has been reduced by 1.1 𝑒𝑉. Moreover, the currents have a positive sign for applied biases higher than −2 𝑉 , indicating that the device is operating under the thermionic conversion regime. Figure B-4. The current voltage characteristics of the in-situ intercalated MWCNT forest under LITEC conditions. The incident laser power is 35 𝑚𝑊. The values of the workfunction of the emitter, the inter-electrode distance and the emitter temperature were calculated to be 3.5 𝑒𝑉, 980 𝜇𝑚 and 𝑇𝐸 = 1530 𝐾, respectively. 123 Appendix C: Publications that are not directly related to the results of this thesis I have helped and contributed to the following Journal paper and conference presentations during my PhD studies: Journal paper: 1. A. T. Koch, A. H. Khoshaman, H. D. E. Fan, G. A. Sawatzky, and A. Nojeh, “Graphenylene Nanotubes,“ The Journal of Physical Chemistry Letters, vol. 6, no. 19, pp. 3982–3987, Oct. 2015 Conference papers: 1. A.H. Khoshaman, A. Koch, T. Pan, A. Nojeh, ”Controlled deposition of SWNT film suspended over nanostructured dielectric,“ The 16th International Conference on the Science and Application of Nanotubes (NT15) 2. A.H. Khoshaman, A. Koch, T. Pan, A. Nojeh, ”Extremely High Electron Transparency of Carbon Nanotube films,“ The 16th International Conference on the Science and Application of Nanotubes (NT15) 3. A.H. Khoshaman and A. Nojeh, ”Modeling of Thermionic Converters through Self-consistent Solution of Vlasov and Poisson equations,“ 28th International Vacuum Nanoelectronics Conference (IVNC 2015) 4. A.H. Khoshaman, A. Koch, A. Nojeh, ”Quantification of Carbon Nanotube Film Properties from Scanning Electron Microscope Images,“ The 16th International Conference on the Science and Application of Nanotubes (NT15) 5. K. Dridi, A. H. Khoshaman, A. Nojeh, G. A. Sawatzky, “Towards compact solar thermionic converters based on carbon nanotubes forests,” The International Vacuum Nanoelectronics Conference (IVNC), Vancouver. 6. H. Fan, A.H. Khoshaman, A.T. Koch, A. Nojeh, “Thermoelectric devices based on sandwiched carbon nanotube forests and the effect of contact electrodes,“ NT14: The Fifteenth International Conference on the Science and Application of Nanotubes, Los Angeles, California, June 2014. 7. H. Fan, A.H. Khoshaman, A.T. Koch, A. Nojeh, “Thermoelectric Performance of Vertically-aligned Carbon Nanotube Forests,“ 97th Canadian Chemistry Conference and Exhibition, Vancouver, British Columbia, June 2014. 8. M. Chang, M.V. Moghaddam, A.H. Khoshaman and A. Nojeh, ”High Temperature Gradient in a Conductor: Carbon Nanotube Forest under the “Heat Trap” Condition,“ EPIBN 2013 124 9. A. Koch, S. Motavas, A. Khoshaman, H. Fan, G. Sawatzky, A. Nojeh, ”Graphenylene-based nanotubes,“ Contributed talk and poster, 15th International Conference on the Science and Application of Nanotubes (NT14), Los Angeles, CA, June 2014.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Self-consistent Vlasov-Poisson analysis of carrier...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Self-consistent Vlasov-Poisson analysis of carrier transport in vacuum-based thermionic/thermoelectronic… Khoshaman, Amir Hossein 2016
pdf
Page Metadata
Item Metadata
Title | Self-consistent Vlasov-Poisson analysis of carrier transport in vacuum-based thermionic/thermoelectronic devices |
Creator |
Khoshaman, Amir Hossein |
Publisher | University of British Columbia |
Date Issued | 2016 |
Description | Thermionic conversion involves the direct conversion of heat, including light-induced heat, from a source such as solar energy, to electricity. The progress of thermionic converters has been limited by issues such as the space charge effect and lack of materials with desirable mechanical, electrical and thermal properties. Nanotechnology could help address some of the main challenges that thermionic converters encounter. However, existing models, which were developed for macroscopic converters, are not adequate for many aspects of nanostructured devices. The work presented in this thesis primarily advances a new model to partially address this void and study emergent thermionic devices. We demonstrate a self-consistent and iterative approach to the Vlasov-Poisson system that overcomes the inherent limitations of the traditional methods. This approach serves as the foundation for more advanced and yet crucial cases of the operation of thermionic converters in the presence of back-emission, grids in the inter-electrode region and low-pressure plasmas. We develop the physics of the device in the presence of grids and demonstrate that momentum gaps could arise in the phase space of the electrons; taking into account these gaps, which had not been noticed in the past, is key to designing efficient thermionic converters and we predict improvements of 3 orders of magnitude in current density using a properly designed grid. We also develop the physics of the device in the presence of low-pressure plasmas, which are prime candidates for reducing space charge. We show that the output power density of a thermionic converter can improve by a factor of ~ 10 using a modest plasma pressure of 500 Pa. On a different front, we have also improved the traditional analytical model and developed an approach to extract the internal device parameters such as emission area and workfunction based on a limited set of experimental output characteristics. These parameters are highly dependent on the operating conditions and ex-situ measurements are not applicable. Therefore, our approach allows for a more systematic study of the device and material properties, which is key to further the development of thermionic converters, in particular based on novel materials and nanostructures. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2016-08-24 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NonCommercial-NoDerivatives 4.0 International |
DOI | 10.14288/1.0308734 |
URI | http://hdl.handle.net/2429/58943 |
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 |
GraduationDate | 2016-09 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nc-nd/4.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2016_september_khoshaman_amir.pdf [ 6.71MB ]
- Metadata
- JSON: 24-1.0308734.json
- JSON-LD: 24-1.0308734-ld.json
- RDF/XML (Pretty): 24-1.0308734-rdf.xml
- RDF/JSON: 24-1.0308734-rdf.json
- Turtle: 24-1.0308734-turtle.txt
- N-Triples: 24-1.0308734-rdf-ntriples.txt
- Original Record: 24-1.0308734-source.json
- Full Text
- 24-1.0308734-fulltext.txt
- Citation
- 24-1.0308734.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}]}"
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.24.1-0308734/manifest