Modelling of a plasma break-down for the SMART tokamak

Master’s Degree Thesis Inter-University Master’s Degree in Nuclear Physics Faculty of Physics, University of Seville Modelling of a plasma break-down for the SMART tokamak Author: Daniel L´opez Aires1 Supervisors: Dr. Manuel Garc´ıa Mu˜noz2 Dr. Carlos Soria del Hoyo3 September 7, 2020 1E-mail addres: [email protected] 2Department of Atomic, Molecular and Nuclear Physics, Faculty of Physics, University of Seville 3Department of Electronics and Electromagnetism, Faculty of Physics, University of Seville Abstract Achieving controlled nuclear fusion on Earth could be a decisive step on the quest towards green energy production since it would bring a virtually renewable energy source without CO2 emissions. The Plasma Science and Fusion Technology group of the University of Seville is currently designing a magnetic fusion device, a spherical tokamak for controlled nuclear fusion research, called SMall Aspect Ratio Tokamak (SMART). This thesis will model the first phase of the initiation of a tokamak (tokamak start-up), the break-down phase, in which the fuel will transition from the gas state to the plasma state, making use of the Fiesta toolbox. The fundamental of tokamak physics and tokamak startup will also be reviewed. The current waveforms of the SMART coilset have been optimized to achieve the desired plasma equilibrium and allow the break-down of the pre-fill gas by reducing the stray poloidal magnetic field and maximizing the loop voltage induced by the inductor solenoid. Several criteria have also been applied to test the feasibility of the breakdown phase, such as the Paschen’s break-down curve, the estimation of the avalanche time, and the calculation of the connection length. The electric potential gained by the electrons as they follow the magnetic field lines has been computed to estimate where the gas will break-down. Break-down of the gas without the use of any auxiliary heating method has been achieved, lasting few milliseconds for gas pressures about 10−4Torr. i Acknowledgements I would like to thank J. Lister, who referred me to Izaskun Garrido, who gave me the main reference to understand the RZIp code. I would also like to thank Geoff Cunningham for his help and comments which help me focus on the right issues regarding break-down, and which also provide us with the Fiesta repository on git and a brief Fiesta documentation. Special thanks to my tutors, this months of work have been really valuable, I have learned a lot, especially about team work in a multidisciplinary environment. I also have to apologize to Carlos for not telling him all the rapid upgrades the SMART reactor was suffering to ensure proper break-down conditions, isolating him a bit from the group and from my work. Thanks to Eli, for her help supervising this document. Special thanks to Scott, who, despite arriving to the group in march, have learned very rapidly and have played a critical role in the changes in the SMART reactor, and to Alessio and Manu too, who also played a critical role in the changes from the engineering point of view. I would also like to thank Jes´us Poley, our talks about fusion and science in general were really interesting and helpful for both of us. Thanks to Antonio for his help with Inkscape. Special mention to the CS-GO team for all the good Friday nights during the Covid-19 quarantine. ii Contents Contents iii 1 Introduction 1 1.1 Nuclear fusion ..................................... 1 1.1.1 Definition of a plasma ............................ 5 1.2 Confinement of charged particles .......................... 6 1.3 Tokamaks ....................................... 9 1.4 Motivation: SMART ................................. 12 1.5 Objectives and outline of the thesis ......................... 13 2 Fundamentals of tokamak physics 15 2.1 Magnetohydrodynamic model of a plasma ..................... 15 2.2 G-S equation ..................................... 16 2.3 Parameters ...................................... 18 2.3.1 Plasma shape and control in tokamaks ................... 18 2.3.2 Safety factor and normalized pressure .................... 20 3 Tokamak start-up 22 3.1 Plasma break-down .................................. 24 3.1.1 Avalanche or break-down time ........................ 28 3.1.2 Voltage and electric field induced by the inductor solenoid ........ 29 4 Fiesta toolbox 31 4.1 Grad-Shafranov solver, EFIT ............................ 31 4.2 RZIp model ...................................... 32 5 Simulation procedure 37 5.1 Field line tracer .................................... 39 iii CONTENTS CONTENTS 5.1.1 Potential .................................... 40 6 Results 41 6.1 Static and dynamic behaviour of SMART ...................... 41 6.1.1 Target equilibrium .............................. 41 6.1.2 Current waveforms .............................. 42 6.1.3 Dynamic simulations ............................. 45 6.2 Break-down results .................................. 46 6.2.1 Magnetic fields ................................ 46 6.2.2 Connection length .............................. 48 6.2.3 Paschen curve and avalanche time ...................... 53 6.3 Discussion ....................................... 57 7 Summary and conclusions 62 8 Bibliography 63 A Calculation of the voltage induced by a solenoid of finite width 68 B State space representation 70 iv Chapter 1 Introduction 1.1 Nuclear fusion as an energy source Nowadays, the human kind is beginning to understand the damage its activity is causing on Earth, that could lead to the destruction of the planet we live in and, as a consequence, of ourselves. A radical change is needed in human’s life before it is too late. One fundamental step is to stop using fossil fuels as an energy source, and use renewable sources instead, like wind energy, solar energy, or geothermal energy. However, there is another energy source, virtually unlimited that could provide a huge step forward this transition, nuclear fusion. Nuclear fusion is a type of nuclear reaction in which two or more atomic nuclei (reactants) Xand Yinteract and produce a heavier nuclei, generally in an excited state (X+Y)∗. This compound nuclei could de-excitate by emiting electromagnetic radiation and, if the excitation energy is sufficiently high, by releasing neutrons (evaporation). X+Y→(X+Y)∗(1.1) Applying the conservation of energy to the reaction, using the laboratory reference frame EX+EY=E(X+Y)∗=E(X+Y)+Eexc ⇒Ti+ (mX+mY)c2=Tf+ (mX+mY)c2+Eexc ⇒Tf−Ti≡Q=−Eexc, (1.2) Since the excitation energy Eexc is always positive, the Qfactor of the reaction is negative, meaning that there is a reference frame in which Tf= 0, but Tican never be zero. This means that this reaction is an endothermic reaction, it needs energy to take place. 1 1.1. NUCLEAR FUSION CHAPTER 1. INTRODUCTION Figure 1.1. Coulomb barrier between two nuclei of mass number Aand a.e=e/(4πε0). Source: [1]. The need of energy to produce a nuclear fusion reaction can be easily understood. Due to the positive charges of the nuclei, their Coulomb interaction is repulsive. However, not always the interaction is repulsive, if the nuclei are close enough, the nuclear interaction appears, and since its much more intense than the electromagnetic interaction, the dominant interaction is nuclear, which will attracts the nuclei, enabling them to approach enough so they can fuse into a new nuclei. If we plotted the potential between two nuclei, it would look similar to the one on figure 1.1. There are two regions, the region for r >>, in which the nuclei are far away from each other so there is no nuclear interaction between them and the potential is the Coulomb potential, and the r << region, in which the nuclear interaction appears, so the potential is the nuclear potential and the nuclei attract each other. The point where Coulomb interaction is compensated by the nuclear interaction is usually estimated as RN≃1.45(A1/3 1+A1/3 2)fm, where A1,A2are the mass number of the nuclei (Aand ain figure 1.1). The main challenge of nuclear fusion reactions is that the Coulomb barrier have to be overcomed, allowing the nuclei to approach enough so they can fuse into a new nuclei. Although nuclear fusion reactions need energy to take place, the reactions could release more energy than the needed to stimulate it. This can be understood taking into account the concept of binding energy B(N, Z), which is the energy needed to split a nuclei into its components, B(A, Z)≡Zmp+ (A−Z)mn−[M(A, Z)−Zme])c2,(1.3) where M(A, Z) is the mass of an atom of atomic number Zand mass number A, and mp,mn and meare the proton, neutron and electron masses respectively. If B/A is plotted, figure 1.2 is obtained. B/A increases with A up to 56Fe, and then it starts decreasing. This means that if combining elements to the left of 56Fe, the compound nuclei is more stable than the initial nuclei, and the difference of binding energy from the stable compound nuclei and the less stable separate nuclei is released in the form of kinetic energy. Equivalently, if a nuclei heavier than 56Fe splits into lighter nuclei, since the final nuclei will have higher binding energy than the original nuclei, the difference of binding energy will be released. 2 1.1. NUCLEAR FUSION CHAPTER 1. INTRODUCTION Figure 1.2. Binding energy per nucleon. The maximum B/A correspond to 56Fe, the most stable element. Source: [2]. Figure 1.2 essentially explains why energy can be obtained from fission reactions, which have lead to the creation of nuclear power plants to obtain energy by fission reactions. But, in the same way it suggest that we could obtain energy as well by pursuing nuclear fusion reaction. In both cases, the fundamental condition to obtain energy is that the energy released is greater than the energy applied. In energy context, it is defined a variable called Qfactor which is the ratio between the power obtained and the power applied to the system, Q=Pobt Papp .(1.4) How can we achieve nuclear fusion reactions on Earth so it can be used as an energy source (controlled nuclear fusion)? The first challenge is that for nuclear fusion to happen, the Coulomb barrier needs to be overcomed (actually, considering quantum tunneling, energy lower than the needed to surpass the Coulomb barrier would be needed to allow fusion reactions). Nuclear fusion reactions considered to be pursued on Earth involve Deuterium, 2 1H because it is an abundant element, it exist on Earth’s oceans comprising 0.015 atom percent of the hydrogen in sea water with the volume of about 1.35 ·109km3[3] (section 1.3). The possible reactions are: 1. 2 1H + 2 1H→3 1H(1.01MeV) + p(3.03MeV), 2. 2 1H + 2 1H→3 2He(0.82MeV) + n(2.45MeV), 3. 2 1H + 3 1H→4 2He(3.52MeV) + n(14.06MeV), 4. 2 1H + 3 2He →4 2He(3.67MeV) + p(14.67MeV), where the kinetic energy each product carries is also indicated. Their cross-sections are plotted in figure 1.3, where the cross-section of the two 2 1H-2 1H reactions have been added, and the X-axis is the projectile energy, 2 1H, assuming the target nuclei at rest. 2 1H + 3 1H reaction is the best option since its cross-section is the highest, and it peaks at the lowest energy1. 1Tritium, 3 1H, do not exist naturally on Earth’s, but it could be produced by certain nuclear reactions with 3 1.3. TOKAMAKS CHAPTER 1. INTRODUCTION a∇Bdrift, that drift the ions downward and the electrons upward, since (1.14) depends on the charge q. This charge separation will create an electric field perpendicular to the magnetic field, so the particles will experience a ~ E∧~ Bdrift, that would drift outward both ions and electrons, provided that this drift does not depend on the charge. These drifts are shown in figure 1.7 (a). (a) Drifts in a torus. Source: [5]. (b) Toroidal (blue) and poloidal (red) directions of a torus. Source: google images, 2019. Figure 1.7. Drifts in a torus and definition of the poloidal and toroidal directions on a torus. To overcome the drifts discussed above, Tokamaks confine the particles by twisting the magnetic field lines. For doing that, in addition to the toroidal field of the torus, a poloidal field is added (field in the poloidal direction, see figure 1.7 (b)). This poloidal magnetic field is created by the plasma itself. Tokamaks need additional coils for controlling the plasma. The plasma itself tend to move radially outward due to poloidal field created by the plasma, which is greater in the inboard region than in the outward, and due to the toroidal shape of the plasma, the plasma pressure also creates an outward force. This outward force is called hoop force (see [11] for an illustrative explanation, and [3] for a rigurous treatment). This coils are called poloidal magnetic field coils (PF coils) since its role is to create poloidal field whose ~ J∧~ Bforce balance the hoop force. For doing that, the current flowing in this PF coils need to flow in the opposite direction to the plasma current. In addition, this coils also help create the elongated shape of tokamak plasmas. An additional set of coils called divertor coils are often used to created a diverted shape in the plasma, which will be explained later. Figure 1.8 shows a sketch of a tokamak with its basics elements. The plasma is contained in the vacuum vessel (VV) (grey coloured in the figure). The toroidal field (green arrows) is created by the toroidal magnetic field coils, and the poloidal field is created mainly by the plasma itself, 10 1.3. TOKAMAKS CHAPTER 1. INTRODUCTION and by the PF coils (plasma shaping). In the center of the device there is a transformer coil that induces a toroidal current in the plasma (red arrows) that creates the poloidal magnetic field. Figure 1.8. Sketch of a tokamak, showing its basic elements, the field lines and the plasma current. Source: google images, 2019. The resulting field lines are helical lines (yellow arrows), which confine most of the plasma particles. Since the plasma current is inductive, tokamaks operate in a pulsed regime. One the challenges of this device is the control of the plasma, and the electromagnetic instabilities that could arise in the plasma, which could lead to the loss of the plasma energy. This events are called Disruptions. For the initialization of a tokamak discharge, called tokamak start-up, the toroidal magnetic field must be previously stablished and the VV is filled with a gas. After the inductor coil induced the electric field by changing its current, the gas will be ionized creating a plasma (this is called break-down). This plasma will start to create the poloidal field that confines the particles. In the meantime, the plasma needs to be heated and the PF coils will be controlling its shape. There are several methods of heating the plasma; to begin with, the plasma current will heat the plasma due to Joule’s effect, which is called ohmic heating. External methods of heating could be the use of electromagnetic waves (the electromagnetic waves will create oscillations of the plasma particles, increasing their energy), injection of neutral particles (the injected particles will accelerate the plasma particles by collisions with them, increasing their temperature), and many more (see chapter 5 of [4] for a description of heating methods6). Finally, the plasma will arrive at the desired configuration with the desired plasma current and shape, ending the start-up phase. 6The websites of currently operating tokamaks also provide information about the way they heat their plasmas. 11 1.4. MOTIVATION: SMART CHAPTER 1. INTRODUCTION 1.4 Motivation: SMART, a SMall Aspect Ratio Tokamak for the University of Seville The Plasma Physics and Fusion Technology Group of the University of Seville is designing a tokamak that will be operating the next year7. Its name will be SMall Aspect Ratio Tokamak (SMART), which reveal the main characteristic of the device, it will be a spherical tokamak rather than a standard tokamak. This reactor will not make fusion reactions; instead, as well as all the existing small and medium size reactors, it will focus on doing tokamak physics research such as plasma confinement, shape, instabilities, etc. This information will be used by large international facilities such as JET or ITER, which will make fusion reactions. The SMART missions are  Study plasma transport and confinement in positive and negative triangularities.  Develop novel diagnostic and control schemes.  Examine electromagnetic stability and control of energetic particles.  Train next generation of fusion physicists and engineers. The main difference between a spherical tokamak and a tokamak is the aspect ratio of the device, which is the ratio of the major radius and the minor radius of the device. If the aspect ratio of the tokamak is <2, the device is called spherical tokamak. Figure 1.9 shows a standard tokamak and a spherical tokamak. Figure 1.9. Tokamaks and spherical tokamaks. Source: [12]. Spherical tokamaks are a desirable approach to controlled nuclear fusion because they are more compact than regular tokamaks, which means lower costs, and spherical tokamak’s plasmas displays better plasma properties such as the so called safety factor qand the β, which will be explained later (reviews of spherical tokamaks features vs regular tokamaks features can be found on [12,13]). Three operational phases have been designed for the SMART tokamak: 7http://www.psft.eu/. 12 1.5. OBJECTIVES AND OUTLINE OF THE THESIS CHAPTER 1. INTRODUCTION  Phase 1. First plasma and proof-of-concept.  Phase 2. Inclusion of Neutral Beam Injection (NBI) heating system. The goal of this phase is the demonstration of plasma shaping.  Phase 3. This phase will explore fusion-relevant operations. The main parameters of this three phases are shown in table 1.1. This parameters will be explained in the following sections. A 3D model of the SMART tokamak is shown on figure 1.10. The coilset configuration will be the same for the three phases. (a) Reactor. (b) Coilset. Figure 1.10. 3D plots of the SMART reactor. The coilset is symmetric with respect to the Z= 0 plane, and several coils are inside the vacuum vessel (VV). The ports in the VV are for plasma diagnosis and for plasma heating methods such as NBI and ECRH, which will be explained later. Source: [14]. 1.5 Objectives and outline of the thesis The goal of this work is to model the initial phase of the tokamak start-up, the break-down of the gas, of SMART, optimizing the device so that the gas breaks-down and turns into a plasma for the first operational phase and the first upgrade, phase two, without any additional heating method. This thesis is organized as follows: chapter 2 presents the basics of tokamak physics, the simplest model of the plasma, and the fundamental equation to describe the tokamak equilibrium, the Grad-Shafranov equation, as well as the characterization of a tokamak plasma. Chapter 3 reviews the tokamak start-up, focusing on the break-down phase, introducing the 13 1.5. OBJECTIVES AND OUTLINE OF THE THESIS CHAPTER 1. INTRODUCTION SMART Phase 1 Phase 2 Phase 3 VV radius(m) 0.8 VV height(m) 1.6 Major plasma radius(m) 0.42 Minor plasma radius(m) 0.24 Plasma elongation κ κ < 2.1 Plasma triangularity δ−0.51 < δ < 0.44 Plasma current(kA) 35 100 500 Toroidal field(T) 0.1 0.3 1.0 Flat-top time(ms) 20 100 500 External heating (KW) ECRH 6 [2.4GHz] 6 [7.5GHz] 200 [- GHz] NBI - 600 600 Table 1.1. Parameters of the phases of the SMART tokamak. Toroidal field is the toroidal field value at the plasma magnetic axis. basic physics of this phase, as well as several criteria to ensure a proper break-down of the gas. Chapter 4 reviews the Fiesta toolbox used for computing the plasma equilibrium and the dynamic behaviour of the plasma using the RZIp model. Chapter 5 summarizes the simulation procedure followed in this work. Chapter 6 contains the simulation results, as well as the discussion of the results with other operating tokamaks. Finally, chapter 7 summarized the work carried out in this thesis, and discuss future work. 14 Chapter 2 Fundamentals of tokamak physics In this chapter the fundamentals of tokamak physics will be reviewed. The first section includes the plasma model used in tokamak physics, the second section includes the derivation of the equation of the tokamak equilibrium configuration, and the third section reviews some basic tokamak parameters. 2.1 Magnetohydrodynamic model of a plasma The Magnetohydrodynamic (MHD) model is a single fluid description of a plasma, i.e., it is a model that treats a plasma like a continuum matter, rather than as a set of particles. This is one of the simplest models to study a plasma, and it assumes several hypothesis, like quasineutrality or negligible electron inertia (see any book about plasma physics for further details, like [5]). The set of equations are ∂ρ ∂t +∇ · (ρ~v) = 0,[Mass conservation] (2.1) ρh∂~v ∂t + (~v · ∇)~vi=~ j∧~ B− ∇p, [Momentum conservation] (2.2) ~ E+~v ∧~ B=η~ j, [Ohm’s law] (2.3) ∂ ∂t(pρ−γ)+(~v · ∇)(pρ−γ) = 0,[Adiabatic behaviour] (2.4) ∇ ∧ ~ E=−∂~ B ∂t ,(2.5) ∇ ∧ ~ B=µ0~ j, (2.6) ∇ · ~ B= 0,(2.7) where ρis the plasma density, ~v its velocity, ηits resistivity, pits pressure (in general the pressure is a tensor, but for this simplified model, it is considered an scalar magnitude), ~ j 15 2.2. G-S EQUATION CHAPTER 2. FUNDAMENTALS OF TOKAMAK PHYSICS its current density, γis the adiabatic index, and ~ Eand ~ Bthe electric and magnetic field the generated by the plasma. Note that the last three equations are a quasi-static limit of Maxwell’s equations. 2.2 Grad-Shafranov equation Figure 2.1(a) displays the coordinate system of toroidal devices. In an equilibrium situation, the magnetic field of a Tokamak produces an infinite set of nested toroidal magnetic flux surfaces1 as shown in figure 2.2, and the magnetic field lines follow an helical path on them as they wind round the torus. The poloidal flux Ψ and the toroidal flux Φ between two magnetic surfaces are defined by dΨ≡~ B·~ dSθ, dΦ≡~ B·~ dSφ,(2.8) where dSθ,dSφare the poloidal and toroidal surface elements, whose magnitude are defined in figure 2.1(b) (its unitary vector is perpendicular to the surface, and the sign is arbitrary, as usually in the magnetic fluxes), and ~ Bis the magnetic field. The basic condition for the (a) Cylindrical coordinate system used in devices with torodial symmetry, (R, φ, Z). R0is called the major radius of the torus, ris called the minor radius. The circumference R=R0defines the toroidal or magnetic axis. Source: http://fusionwiki.ciemat.es/wiki/ Toroidal_coordinates. (b) Toroidal (T) and poloidal (P) surface elements between two magnetic flux surfaces. Source: [15]. Figure 2.1. Cylindrical coordinate system for toroidal devices, and definition of the poloidal flux in a torus. equilibrium is that the force on the plasma be zero at all points, so the momentum conservation 1A given surface is a magnetic flux surface if it satisfies ~ B·~n = 0, where ~n is the normal vector of the surface. That is, the magnetic field do not cross the surface. This is only a visual way to understand the magnetic field, since there would be an infinite number of magnetic flux surfaces inside a tokamak. 16 2.2. G-S EQUATION CHAPTER 2. FUNDAMENTALS OF TOKAMAK PHYSICS equation (2.2) leads to ~ j∧~ B=∇p. (2.9) This implies ~ B· ∇p= 0, so there is no presure gradient along the magnetic field lines, which means the magnetic surfaces are also pressure surfaces. (2.9) also implies ~ j· ∇p= 0, and as a consequence the current lie in the magnetic surfaces. In what follows the Grad-Shafranov equation, one of the most fundamental equations of MHD equilibrium, will be derived. The idea for this equations is to rewrite (2.9) to have a scalar equation instead of a vector equation. From (2.7), taking into account the axysymmetry, and using the coordinate system of figure 2.1(a), setting R0≡0, 1 R ∂(RBR) ∂R +∂BZ ∂Z = 0.(2.10) The function of that scalar equation will be the function ψ, called the stream function, which is defined as ψ≡RAφ, where Aφis the toroidal component of the vector potential ~ A. With this function, the poloidal magnetic field can be written as BR=−1 R ∂ψ ∂Z , BZ=1 R ∂ψ ∂R,            ⇔~ Bθ=1 R∇ψ∧∧ φ, (2.11) where ∧ φis the toroidal unit vector (the magnetic field can be expressed as ~ B=~ Bθ+~ Bφ). It can be shown that Ψ = 2πψ [15] (section 6.2). It is usual to label the magnetic surfaces with ψ, also called the magnetic flux. This means that p=p(ψ), since magnetic surfaces are also pressure surfaces. From the symmetry of ~ j, it can be introduced a function fthat verifies jR=−1 R ∂f ∂Z , jZ=1 R ∂f ∂R,            ⇔~ jθ=1 R∇f∧∧ φ . (2.12) Comparing (2.12) with (2.6) leads to f=RBφ µ0 ,(2.13) where µ0is the vacuum magnetic permeability and the subscript φindicates the toroidal component. It can be shown that fis a function of ψ[4] (section 2.3). Equation (2.9) can be expanded as ~ jθ∧∧ φ Bφ+jφ ∧ φ∧~ Bθ=∇p, (2.14) 17 2.3. PARAMETERS CHAPTER 2. FUNDAMENTALS OF TOKAMAK PHYSICS where jφ,~ jθare the magnitude of the toroidal current, and the poloidal current density vector respectively. Substituting (2.12) and (2.11) into (2.14), we get, using that ∧ φ·∇ψ=∧ φ·∇p= 0 (consequence of the toroidal symmetry) Bφ R∇f+jφ R∇ψ=∇p. (2.15) Now, applying the chain rule on ∇fand ∇p, ∇f=df dψ∇ψ, ∇p=dp dψ∇ψ,      (2.16) Introducing (2.13) and (2.16) into (2.15) lead to: −µ0f R2 df dψ∇ψ+jφ R∇ψ=dp dψ∇ψ⇒jφ=µ0f R df dψ +Rdp dψ.(2.17) ∇ψcan be removed from (2.17) since ∇ψ= 0 correspond to the trivial solution. To get jφas a function of ψ, we substitute (2.11) on (2.6), obtaining µ0~ j=µ0jφ ∧ φ+1 R∇(RBφ)∧∧ φ⇒ −µ0Rjφ=R∂ ∂R1 R ∂ψ ∂R+∂2ψ ∂Z2,(2.18) Finally, if we substitute (2.18) on (2.17), we get the Grad-Shafranov equation, R∂ ∂R1 R ∂ψ ∂R+∂2ψ ∂Z2=−µ0R2dp(ψ) dψ −µ2 0f(ψ)df(ψ) dψ .(2.19) (2.19) is one of the fundamental equations of MHD equilibrium. It is a second order partial differential equation that calculates the equilibrium in toroidal devices, given the functions p(ψ) and f(ψ). In figure 2.2 (b), we can see a typical solution of this equation, showing that the surfaces are shifted with respect to the magnetic axis, which is the major radius of the innermost surface. 2.3 Tokamaks parameters The most important parameters of a tokamak equilibrium, as well as about plasma shape and control will be introduced in this section. 2.3.1 Plasma shape and control in tokamaks The first concept that must be introduced is the plasma boundary. The boundary of the plasma is the outermost closed magnetic surface contained in the VV, called Last Closed magnetic Flux 18 2.3. PARAMETERS CHAPTER 2. FUNDAMENTALS OF TOKAMAK PHYSICS (a) Magnetic flux surfaces of a tokamak equilibrium, forming a set of nested cylindrical surfaces. Source: [4]. (b) Typical solution of the Grad- Shafranov equation. Source: [4]. Figure 2.2. Magnetic flux surfaces of a tokamak equilibrium, and typical solution of the Grad-Shafranov equation. Surface, LCFS. The particles inside this outermost surface follow the field lines that remain in the plasma, but particles that follow the external field lines will end up escaping from the plasma and colliding with the VV (the particles follow the magnetic field lines, but if the magnetic field lines are not closed inside the VV, particles will collide with the VV as they follow them). The boundary can be created by a set of coils or by a limiting material that will touch the plasma (either the vacuum vessel or another element). The first method is to create the LCFS displaying one or more X-points by using a set of coils. X-points are saddle points where ∂ψ ∂Z =∂ψ ∂R = 0, so the poloidal magnetic field is zero (see (2.11)). The outermost closed surface is called separatrix, and a plasma confined this way is called a diverted plasma, and the coils used to create it are called divertor coils. The second method is to limit the plasma by the VV or an specific material, so that the plasma is touching that material. The plasma is then called a limited plasma. In figure 2.3 (a) a limited plasma and a diverted plasma with one X-point is showed. For shape control of the plasma, the following parameters are introduced to describe the 19 3.1. PLASMA BREAK-DOWN CHAPTER 3. TOKAMAK START-UP with each other. The connection length can be calculated numerically by integrating the magnetic field lines equation, but also an empirical formula is used to estimate it in the poloidal field null region [25]: L≃0.25aeff Bϕ(Rnull) < Bθ>,(3.5) where Bϕ(Rnull) is the toroidal magnetic field at the center of this region Rnull and < Bθ> is the average poloidal field in the surface of this region. Regarding aeff, there are two visions in the bibliography: i) aeff is the linear distance to the closest wall. Sometimes this is estimated as the minor radius of the plasma target equilibrium configuration [21,23]. ii) aeff is the minor radius of the field null region [20,25,19]. The variation in the electron density considering both ionization and losses is dne dt =ne(νion −νloss).(3.6) The ionization rate can be assumed constant since electrons will acquire a constant speed quickly after the start of the ionization process. About the loss rate, the connection length may vary since as the toroidal electric field is induced inside the VV, eddy currents are also induced in the VV, which will affect the field, and hence alter the connection length L. However, as a first approximation, we could assume Lto be constant. In that case the electron density as a function of time will be ne(t) = ne(0) exp [(νion −νloss)t],(3.7) where ne(0) is the electron density at t= 0, the time at which the toroidal electric field is induced. For a proper ionization of the gas, the ionization rate must be greater than the loss rate, so that the electron density increases. Imposing this, and introducing (3.1) and (3.4) in (3.6) gives dne dt =nevk(α−1/L)>0⇒(α−1/L)>0⇒αL > 1,(3.8) taking into account vkis the module of the speed, ergo, positive. For Townsend avalanche to proceed, αL > 1 is needed. Setting αL = 1 will give the limit condition so that the avalanche nor increase nor decrease. Introducing this condition on (3.2) gives Eϕmin =C2p ln(C1pL),(3.9) where Eϕmin is the electric field needed for this condition. Since this condition do not ensures avalanche, [25] states that Eϕ>2Eϕmin for a reliable start-up, so that the avalanche increases. 26 3.1. PLASMA BREAK-DOWN CHAPTER 3. TOKAMAK START-UP 10 -5 10 -4 10 -3 p (Torr) 10 -1 10 0 10 1 10 2 10 3 10 4 Emin (V/m) Paschen's curve for Hydrogen L=10m L=50m L=100m Figure 3.2. Paschen’s curve for for H2as pre-fill gas, showing different connection lengths. The plot of Eϕmin as a function of pfor a given Lis called the Paschen’s curve. An example is shown on figure 3.2, for H2as pre-fill gas. The constant are C1= 510m−1Torr−1 and C2= 1.25 ·104Vm−1Torr−1. It can be seen from the figure that all the lines displays a minimum electric field needed for a certain pressure. This can be easily understood: for high pressures, the mean free-path will be too short (the mean free path λis given by λ= 1/(C1p) [21]), so for the electrons to gain enough energy to ionize, the electric field needs to be high. For low pressures, the mean free-path will be too long meaning that there would be too little collisions before the electrons are lost, although all the collisions will be ionizing collisions because the electron will have gained enough energy. This results in a vertical asymptote, so that to the left of that asymptote the loss rate is greater than the ionization rate, regardless of the field, because there is not enough molecules to create avalanches. A widely used empirical criteria for a reliable startup is the so-called Lloyd’s criteria [24] Eϕ Bϕ Bθ >1000Vm−1.(3.10) In the case of plasma break-down assisted by Electron Ciclotron Resonance Heating (ECRH), which is a heating method based on the irradiation of electromagnetic waves of certain frequencies to the ionized gas so that it absorbs them helping ionize the gas [19], the minimum value is 100Vm−1instead of 1000Vm−1[26]. In the DIII-D tokamak3it has been shown [27] that the break-down of the pre-fill gas did not occur where the stray poloidal field was minimum; instead, it took place where the potential gained by the electrons along their path Rfield line ~ E·~ dl was greatest. This potential is not the electrostatic potential because the electric field has been created by a changing magnetic field. An electrostatic field ~ Ecan be written as ~ E=−∇V, where Vis the electrostatic potential. The line integral between two points Aand Bis RB A~ E·~ dl =−V(B) + V(A), so it does not depend on the path followed, only on the starting and ending points. In the presence of transient magnetic fields, the electric field is given by ~ E=−∇V−∂~ A/∂t, where ~ Ais the vector 3https://www.ga.com/magnetic-fusion/diii-d. 27 3.1. PLASMA BREAK-DOWN CHAPTER 3. TOKAMAK START-UP potential. The line integral in this case would be RB A~ E·~ dl =RB A−∇V·~ dl −RB A∂~ A/∂t ·~ dl= −V(B) + V(A)−RB A∂~ A/∂t ·~ dl. To compute the last integral, the path followed from Ato B needs to be known, and hence, the final result will depend on the path followed. This means that an electron gains energy after doing a revolution inside the tokamak returning to the starting point. The NSTX tokamak4also confirms this empirical results [28]. This seems understandable since to create avalanches, the electrons have to gained enough energy to ionize, so it is not only relevant that they follow a path with the minimum deviation to collide as much as possible, but that they gain as much energy as possible in its path so that most of the collisions are ionizing collisions. 3.1.1 Avalanche or break-down time When Coulomb collisions dominate over neutral-electron collisions, the fuel starts to behave like a plasma, entering the burn-through phase. This phase is reached usually when the ionization fraction of the gas is about 5% [19]. At this stage, the avalanche as described above, with a gas being ionized by a toroidal electric field stops being valid, and the further ionization is made by the plasma itself, if the power losses are counterbalanced. The time when the this phase starts (or when the break-down phase ends) can be estimated using eq (3.7). Introducing the concept of ionization fraction fi, which is the ratio between the electron density created by ionization and the neutral pre-fill gas density npre, which is fi≡ne/2 npre ,(3.11) where the 2 accounts for the 2 electrons that H2(and also He) gives. The pre-fill gas density is related to the gas pressure by the ideal gas law, p=npreKBTpre. Substituting this into (3.7) yields, assuming ne(0)=1 which is an standard assumption [19,23], fi(t) = 1 2npre exp [(νion −νloss)t] = 1 2 KBTpre pexp hvk(α−1 L)ti.(3.12) Setting fi= 5% on (3.12) will give an estimation for the time needed for the avalanche or break-down phase to end, the avalanche time tava: tava = ln 2·0.05 ·p KBTpre  vk(α−1 L) = ln 2·0.05 ·p KBTpre  43Eϕ phC1pexp −C2p Eϕ−1 Li.(3.13) Note that it is positive since α−1/L > 0 which is the condition for the avalanche to occur. 4https://www.pppl.gov/nstx. 28 3.1. PLASMA BREAK-DOWN CHAPTER 3. TOKAMAK START-UP 3.1.2 Voltage and electric field induced by the inductor solenoid To start the break-down phase, the inductor solenoid (Sol) is pre-charged with a certain current, and then its current its ramped down rapidly to induce the electric field. In this subsection the electric field and the voltage induced will be derived. The voltage induced by the ramp down of the Sol, called loop voltage, can be easily computed using Faraday’s law in its integral form: ε≡Vloop =−d dt ZS ~ B·~ndS, (3.14) where Sis the surface enclosed by the loop, which is a circle of an arbitrary radius, and ~n its unit vector, which will go in the Zdirection. The magnetic field ~ Bof a solenoid also goes in the Zdirection, so the loop voltage is, choosing ~n =∧ Z: Vloop =−d dt ZS BSoldS =−ZS dBSol dt dS. (3.15) The derivative can be introduced into the integral since the integration variables do not vary over time. The only thing that vary with time is the Sol current. Its derivative can be computed easily since the Sol current is decreased linearly, so it satisfies ISol(t) = mt +n. (3.16) At t= 0, ISol ≡I0, and the ramp goes down until ISol(t1)≡I1(I0> I1), so the Sol current is ISol(t) = I1−I0 t1 t+I0.(3.17) Solving (3.15) for the magnetic field of a solenoid of inner radius Rin and outer radius Rout > Rin with N/l turns per unit length, ignoring border effects (i.e., considering infinite length), gives (see appendix A) Vloop =−µ0 N l I1−I0 thπR2 in +π(R2 out −R2 in)−2π Rout −Rin 1 3(R3 out −R3 in)−Rin 2(R2 out −R2 in)i. (3.18) Since the slope of the Sol current is negative, the loop voltage induced is positive. Once the loop voltage has been computed, the electric field can also be computed. For computing it, the integral form of Faraday’s law will be used: Vloop =ε≡IΓ ~ E·~ dl =−d dt ZS ~ B·~ndS, (3.19) where Γ is the loop’s perimeter. To compute the electric field ~ Efrom this equation we must know its symmetry. To do so, let’s take into account that (3.19) is very similar to Amp`ere’s law IΓ ~ B·~ dl =µ0ZS ~ J·~ndS =µ0Ienclosed.(3.20) 29 3.1. PLASMA BREAK-DOWN CHAPTER 3. TOKAMAK START-UP It is widely known that the magnetic field of a infinite cilindral conductor carrying a current density ~ J=J∧ Zis ~ B=B∧ ϕ(see any book of classic electromagnetism, like [29]). And, since the roles of ~ Band ~ Jin Amp`ere’s law are the same as the roles of ~ Eand ~ Brespectively in Faraday’s law, provided that ~ B=B∧ Zin Faraday’s law, the electric field must be ~ E=E∧ ϕ. Furthermore, taking into account the axysymmetry and the assumption of no border effects (infinite length), the electric field can only depend on the radial coordinate R, so ~ E=E(R)∧ ϕ. Considering this, the field can be computed, since ~ dl =Rdϕ ∧ ϕ: Vloop =IΓ ~ E·~ dl = 2πRE(R)⇒E(R) = Vloop 2πR .(3.21) The electric field will be more intense near the inner wall of the VV. The same applies to the toroidal magnetic field, which also goes as 1/R. The potential gained in that region will be the greatest, enhancing the ionization processes, and the magnetic field lines will be mostly toroidal lines in that region (since the toroidal field will be greater in that region, the effect of the poloidal field will be reduced as a consequence), reducing the electron losses. This two factors explain why the pre-fill gas breaks-down near the inner side of the VV in general. 30 Chapter 4 Fiesta toolbox The Fiesta toolbox is an object-oriented toolbox programmed in MATLAB by G. Cunningham at the Culham Centre for Fusion Energy [30,31]. It was created for equilibrium calculations such as solving Grad-Shafranov equation and plasma shape control. However, dynamic calculations are also included, using the RZIp model. The reactor geometry is introduced using its object architecture. 4.1 Grad-Shafranov solver, EFIT Fiesta contains what is called free boundary equilibrium solvers. The main characteristic of this type of solvers is that the plasma boundary is not known (on the contrary of fixed boundary problems) and hence will also be part of the solution. The poloidal magnetic flux ψis divided into two components, the created by the plasma and the created by the external current carrying elements, i.e., the coils and the VV, which will carry eddy currents. This external elements will be called structrure. The Grad-Shafranov equation is solved by an iterative method, according to a given tolerance. In particular, we have used the EFIT solver [32], which is a free-boundary equilibrium solver that can obtain the external currents needed to achieve a plasma equilibrium with plasma parameters closer to the desired plasma parameters. The plasma current profile used is defined by the so-called Topeol2 model,      dp dψ =j0 r0 βθ(1 −ψN), fdf dψ =µ0r0j0(1 −βθ)(1 −ψN), (4.1) where r0is the X-point radial coordinate, j0the plasma current and ψNthe normalized flux, 31 4.2. RZIP MODEL CHAPTER 4. FIESTA TOOLBOX defined as ψN≡ψ−ψaxis ψboundary −ψaxis ,(4.2) where ψaxis and ψboundary are the values of the poloidal flux at the center of the magnetic flux surface and at its boundary respectively. The constants r0and j0are adjusted at each step of the iterative process. The field and the flux at the computational grid are calculated using Green’s functions (see [33] for a review of the Green’s function formalism applied to tokamaks). 4.2 RZIp model The rigid current displacement model or RZIp (R,Zand Ip) [34,35] is a model used to describe the dynamic behaviour of a tokamak equilibrium. The plasma is modelled as a rigid conductor so that the plasma can move radially and vertically, but the plasma shape and current distribution remain constant. Hence, the plasma can be identified by its radial and vertical position, Rgeo and Zgeo, and its current, Ip. The RZIp model is a linearised model based on the assumption that small perturbations leads to small changes of the target equilibrium. Its equations are the circuit equation and the radial and vertical force balance equations. This model considers that a tokamak is made of active and passive elements. An element is considered an active element if is subjected to a external excitation such as current or voltage supplies and passive in any other case. The active structure will then be the coilset, and the passive structure the vacuum vessel. There are several ways of calculate the plasma self-inductance. Fiesta’s RZIp uses this approximate formula [35] Lp=µ0Rgeo4πRgeo |< BZp>| µ0Ip−βθ −βθ−1 2,(4.3) where BZpis the BZcreated by the plasma (in our case, the plasma current flows in the positive toroidal direction, resulting in a negative < BZp>). The four equations of the RZIp model are, indicating matrices with bold letters [34] d(LpIp+MpsIs+RGeoEIp) dt +Ipρp=Vp[Ipequation],(4.4) d(LsIs+MspIp) dt +ΩsIs=Vs[Isequation],(4.5) d(mp˙ Rgeo) dt =1 2I2 p ∂Lp ∂Rgeo +Is ∂Msp ∂Rgeo Ip+EI2 p 2[Rgeo equation],(4.6) d(mp˙ Zgeo) dt =Is ∂Msp ∂Zgeo Ip[Zgeo equation],(4.7) 32 4.2. RZIP MODEL CHAPTER 4. FIESTA TOOLBOX where Eis a constant matrix, Mps is the inductance between the plasma and the structure (both active and passive), Isthe current of the structure, Vsand Vpthe voltage of the structure and the plasma respectively and Ωsthe resistance matrix of the structure. Is is common to neglect the plasma mass mp, which is also what Fiesta’s RZIp does. This set of equations are linearized and cast in the state space form to solve it, whose general structure is (see appendix B)    dx dt =Ax +Bu, y=Cx +Du. (4.8) The linearizations of eq (4.4), (4.5), (4.6) and (4.7) setting mp= 0 result in:                      Ls ∂Mps ∂Zgeo 0 ∂Mps ∂Zgeo 0 ∂2Mps ∂Z2 geo 0 Is0 I0 p ∂Mps ∂Rgeo 0 ∂2Mps ∂Zgeo∂Rgeo 0 Is0 I0 p Mps 0 ∂Msp ∂Rgeo 0Mps0 ∂2Mps ∂Rgeo∂Zgeo 0 Is0 I0 p 0 1 2 ∂2Lp ∂R2 geo 0+∂2Mps ∂R2 geo 0 Is0 I0 p ∂Lp ∂Rgeo 0+∂Mps ∂Rgeo 0 Is0 I0 p +µ0 2πS l2βpR0 geo ∂Lp ∂Rgeo 0+∂Mps ∂Rgeo 0 Is0 I0 p +µ0 2πS l2βpR0 geo L0 p+µ0 2πS l2βpR0 geo                       ∗dx dt + +               Ω0 s0 0 0 0 0 0 0 0 0 0 0 0 0 0 Ω0 p               x=         I0 00 00 01              δVs δVp      , (4.9) where 0indicates equilibrium value, so R0 geo indicates the value at the target equilibrium, Sis 33 4.2. RZIP MODEL CHAPTER 4. FIESTA TOOLBOX the plasma cross-section of the equilibrium configuration, and lits perimeter. The state vector xis x=         Is−Is0 (Zgeo −Z0 geo)I0 p (Rgeo −R0 geo)I0 p Ip−I0 p         =         δ(Is) δ(Zgeo)I0 p δ(Rgeo)I0 p δ(Ip)         .(4.10) and the input vector uis u=  δVs δVp .(4.11) However, Fiesta’s RZIp do not uses the voltages as inputs, it uses the coilset currents as inputs, so (4.9) is inverted to have currents as inputs. The system (4.9) has the form Mdx dt +P x =Qu, and comparing it with the state space system (4.8), we conclude A=−M−1P, B=M−1Q. (4.12) The self and mutual inductances on (4.9) are computed with the Green’s functions of the magnetic field by relating them to the field created by the structure and by the plasma. Using (2.11) and the division of the flux into the created by the plasma and the created by the structure: ψ=ψp+ψs, BR=BRp+BRs, BZ=BZp+BZs, BRp=−1 R ∂ψp ∂Z , BRs=−1 R ∂ψs ∂Z , BZp=1 R ∂ψp ∂R , BZs=1 R ∂ψs ∂R . (4.13) The output vector ycontains the state vector variables and diagnostic measurements such as the poloidal field or flux measurements. In our case, we only use poloidal field measures so 34 4.2. RZIP MODEL CHAPTER 4. FIESTA TOOLBOX the output vector is y=Cx ⇒            δ(Is) δ(Zgeo)I0 p δ(Rgeo)I0 p δ(Ip) Bθn            = =             I0 0 0 01 0 0 00 1 0 00 0 1 ∂(~ B·~n)n ∂Is ∂(~ B·~n)n ∂(δ(ZGeo)I0 p) ∂(~ B·~n)n ∂(δ(RGeo)I0 p) ∂(~ B·~n)n ∂Ip                     δ(Is) δ(Zgeo)I0 p δ(Rgeo)I0 p δ(Ip)         , (4.14) where the subscript nin the poloidal field Bθindicates the measure n, and ~n is the normal vector of the loop that measures the poloidal field [34] (section 3.3). The feed-forward matrix Dis zero provided that there is no direct relation between the inputs and the outputs. Using (4.14), the coil currents needed to obtain poloidal field null in the region where the sensors are placed (null currents) can be computed by simply setting to zero the elements of yrelated to the field, and computing the coil currents needed, using only the coil currents in x, and the corresponding elements of the Cmatrix. This is crucial for the break-down phase, since the poloidal field has to be null out as much as possible to ensure the Townsend avalanche. The system (4.9) is solved by a Runge-Kutta adaptive method (function ode45 from MATLAB), obtaining the time evolution of the plasma current, the eddy currents in the vacuum vessel and the coilset and plasma voltages, given the current profile of the coilset. As stated previously, RZIp is a perturbative model on an equilibrium state. This equilibrium configuration have been obtained with the Grad-Shafranov solver of Fiesta, so RZIp is applied after the equilibrium configuration is obtained. Fiesta is a complex toolbox, of thousands of code lines, and it is scarcely documented, so direct code reading is imperative. However, for this thesis is sufficient to understand the basis of the Fiesta toolbox as described above. Figure 4.1 shows a flow diagram of the Fiesta toolbox as it has been used. 35 6.1. STATIC AND DYNAMIC BEHAVIOUR OF SMART CHAPTER 6. RESULTS displays positive triangularities, about 0.2, although plasmas with negative triangularities can also be obtained with the same coilset configuration if the PF2 coils act as divertor coils, so that they take the role of the Div1 coils, and upper and lower bulks will be placed on the PF2 coils instead of on Div1. The parameters q0 and q95 are the safety factor at a surface with ΨN=0 and 0.95 respectively (ΨNis zero for the magnetic axis, and 1 for the plasma boundary or separatrix). Div2 current is zero in the target equilibrium because Div2 is only used in the break-down phase. (a) Phase 1. (b) Phase 2. Figure 6.1. Target equilibrium for both phases, showing the coil currents for the equilibrium configuration and several plasma parameters at the right. 3 bulk region are showed on both plots (upper, lower and central), being the plasma in the central bulk. The colorbar indicates the poloidal flux Ψ. Eddy currents in the VV are not included here. 6.1.2 Current waveforms As stated on chapter 5, the current waveforms were given as an input lacking on the currents during the break-down, which were calculated using RZIp. Previous to the RZIp simulation, the region where the poloidal field will be nulled out (poloidal null region) had to be chosen. The size and location was optimized manually to maximize the empirical connection length (3.5) and reducing the current needed to create the null region. The best shape and size are the one showed on figure 6.2, an square of 0.30x0.30m2 centered at Rnull = 0.31m. The sensors, or magnetic field probes are shown. The region is placed close to the inner wall since its near the inner wall where the gas breaks-down, as explained on 42 6.1. STATIC AND DYNAMIC BEHAVIOUR OF SMART CHAPTER 6. RESULTS Coil I(kA) Phase 1 Phase 2 Sol -0.15 -0.30 PF1 -0.47 -1.49 PF2 -0.07 -0.14 Div1 0.30 1.00 Div2 0.00 0.00 Table 6.1. Coilset currents for the target equilibrium configuration for both phases. Div2 current is zero because it is used only for the break-down phase. Parameter Phase 1 Phase 2 Rgeo(m) 0.42 0.42 A1.82 1.84 κ1.95 2.01 δ0.20 0.24 q0 1.14 1.03 q95 7.23 7.02 βϕ(%) 2.92 3.71 βθ0.73 0.91 βN(%) 1.98 2.56 Table 6.2. Plasma parameters on the target equilibrium for both phases. Zgeo=0 due to the symmetry with respect t the Z= 0 plane. section 3.1. The full current waveforms, after computing the current needed to null out as much as possible the poloidal field in the sensor region with RZIp, is showed on figure 6.3. This figure also includes the current gradient at each time, which will be denoted with ∆. Each current waveform have 7 dots indicating the phases of the plasma discharge, as indicated in figure 6.3 (b):  Point 1, the starting point: all the currents are zero.  Point 2, pre-pulse: the Sol coil is charged, and the PF and Div coils are ramped on to the currents needed to null out the poloidal field created by the Sol (null currents).  Point 3, 1 º ramp-down: at this point, t=0ms, the Sol current is ramped down, inducing the loop voltage that will ionize the pre-fill gas. The PF and Div coils are the null currents to null out the poloidal field so the avalanche can succeed.  Point 4, 2 º ramp-down: at this point, the Sol coil change its decreasing slope to make it less abrupt. The PF and Div coil currents are still the null currents. The break-down of the pre-fill gas will happen between point 3 and 4. The time interval between them was fixed to match the avalanche time.  Point 5, target equilibrium: at this point the target equilibrium is reached. The PF and Div coil currents are the values at the target equilibrium. 43 6.1. STATIC AND DYNAMIC BEHAVIOUR OF SMART CHAPTER 6. RESULTS Figure 6.2. Sensors for the poloidal null region. In that region the poloidal magnetic field will be nulled out (as much as possible) by the PF and Div coils.  Point 6, end target equilibrium: the target equilibrium is sustained from point 5 to point 6, being the time interval between them the flat-top time, which is 20ms on phase 1 and 100ms on phase 2. The Sol coil current vary so that the desired plasma current is maintained, while the PF and Div currents are kept constant  Point 7, termination of the plasma discharge: all the coil currents are turned off, ending the plasma discharge. To model the plasma current, Spitzer’s resistivity have been used, using Zeff = 2 as effective atomic number to account for Carbon impurities because the plasma-facing material in SMART will be Carbon. Due to the low resistivity on phase 1, the Sol current has to increase from point 5 to 6 to prevent further increase of the plasma current. The Sol waveform have 2 decreasing slopes, the first one with a higher slope to induce the electric field and the following with a lower slope to continue increasing the plasma current but at a lower rate. The reason for this two slopes is because of engineering requirements. With the design of the power supplies for phase 1, the first operational phase, the coilset current gradients could not surpass 50A/ms. But this slope was not enough to allow the break-down of the pre-fill gas (the electric field induced was below the paschen’s curve, so break-down was not possible, the loss rate were greater than the ionization rate), so finally it was decided that there would be a first slope to induce the loop voltage that would surpass the operational limits and would be sustained until the break-down phase ends, and another slope within the operational limits once the break-down phase was finished. Note that the rest of the ramp currents for 44 6.1. STATIC AND DYNAMIC BEHAVIOUR OF SMART CHAPTER 6. RESULTS phase 1 are below 50A/ms. The time interval between point 1 and 2 was chosen so that the current ramp is lower than 50A/ms. The second slope was designed so that the plasma current rises up at rate from 1 to 10MA/s, as was found in spherical tokamaks with ohmic startup [27,38,39,40]. The above discussion applies to phase 1 only, but phase 2 was also built with two slopes since this two slopes Sol waveform is less demanding for the power supplies. -40 -20 0 20 40 60 80 time (ms) -0.5 0 0.5 1 1.5 Iinput (kA) Input currents ph1 Sol PF1 PF2 Div1 Div2 (a) Input currents of phase 1. 0 50 100 150 time (ms) -2 -1 0 1 2 3 4 Iinput (kA) Input currents ph2 Sol PF1 PF2 Div1 Div2 7 6 45 1 2 3 (b) Input currents of phase 2, indicating the discrete points of each of the waveforms. -40 -20 0 20 40 60 80 t (ms) -350 -300 -250 -200 -150 -100 -50 0 50 I (A/ms) (Input currents) ph1 Sol PF1 PF2 Div1 Div2 (c) Currents gradients of phase 1. 0 50 100 150 t (ms) -300 -200 -100 0 100 200 I (A/ms) (Input currents) ph2 Sol PF1 PF2 Div1 Div2 (d) Currents gradients of phase 2. Figure 6.3. Coilset currents given as inputs and current gradients for both phases. The desired plasma current is sustained for 20ms in phase 1 and for 100ms in phase 2 (flat-top time). 6.1.3 Dynamic simulations The RZIp outputs, the plasma current and the net eddy current induced in the VV as a function of time are shown on figure 6.4. The coil waveforms have also been included for a better understanding of the outputs. Regarding the plasma current, it can be seen the different growth rate of it due to the 2 slopes of the Sol current, being the first one more abrupt. The plasma is sustained for the desired time, 20ms and 100ms for phase 1 and 2 respectively. At 45 6.2. BREAK-DOWN RESULTS CHAPTER 6. RESULTS the end of the discharge the plasma current is not zero. To null it, cooling down the plasma by allowing air to enter the VV breaking the vacuum would be a method to null out the plasma current. Regarding the eddy currents, they follow the expected behaviour according to the input currents. The highest values are about 30kA, which is about the plasma current in phase 1. The eddy currents at t=0ms are 1.6kA for phase 1 and 2.5kA for phase 2. Figure 6.5 includes the plasma current and its gradient, showing that the growth rates are kept within the desired range, from 1 to 10MA/s (1kA/ms=1MA/s). 6.2 Break-down results 6.2.1 Magnetic fields The magnetic field at the start of the Sol ramp down (dot 3 in figure 6.3 (b) corresponding to t=0ms) have been calculated using the Grad-Shafranov solver. This calculation includes the eddy currents but not the plasma current, not started yet a this time. The Earth’s field was also added to avoid extremely low values of the poloidal magnetic field, as commented on chapter 5. Figure 6.6 includes plots of the toroidal field Bϕ, the poloidal field Bθand the LLoyd’s critera, EϕBϕ/Bθ. Due to the symmetry of the field with respect to the Z=0 plane, only the upper portion of the VV is showed1. Small gaps are appreciable in the VV corners due to rounding issues, without further influence in the simulations. Bϕis the field created by the toroidal field coils, which decays with the radial position as 1/R. The Bθplots displays a more complicated pattern due to the addition of the eddy currents. Without the eddy currents, the lowest Bθis contained within the polidal field null region (dashed in green), but the addition of the eddys create a more complicated pattern. This is a result of the fact that RZIp can not run using a target equilibrium which contains eddy currents, resulting in that the currents for the poloidal field null configuration (null currents) do not take into account the effect of the eddy currents. Hence, the eddy currents have to be added manually, worsening the field to the point that the lowest poloidal field is not in the poloidal field null region. In the central region of the VV, Bθis about 1.1·10−4T in phase 1 and 2.1·10−4T in phase 2. Regarding Lloyd’s empirical criteria, it can be seen that the criteria for ohmic break-down is satisfied (EϕBϕ/Bθ>1000V/m) near the inner wall in both phases, with a larger region in 1Actually, due to the addition of the Earth’s field, the field is not completely symmetric with respect to the Z=0 plane, but the differences are negligible, about 10−7T. 46 6.2. BREAK-DOWN RESULTS CHAPTER 6. RESULTS -40 -20 0 20 40 60 80 time (ms) -0.5 0 0.5 1 1.5 Iinput (kA) Input currents ph1 Sol PF1 PF2 Div1 Div2 (a) Input currents of phase 1. 0 50 100 150 time (ms) -2 -1 0 1 2 3 4 Iinput (kA) Input currents ph2 Sol PF1 PF2 Div1 Div2 7 6 45 1 2 3 (b) Input currents of phase 2, indicating the discrete points of each of the waveforms. -40 -20 0 20 40 60 80 time (ms) 0 5 10 15 20 25 30 35 40 Ip (kA) Plasma current ph1 (c) Plasma current of phase 1. 0 50 100 150 time (ms) 0 20 40 60 80 100 Ip (kA) Plasma current ph2 (d) Plasma current of phase 2. -40 -20 0 20 40 60 80 time (ms) -10 -5 0 5 10 15 20 25 30 IVV (kA) Net eddy current in the VV ph1 (e) Net eddy current in the VV of phase 1. 0 50 100 150 time (ms) -40 -30 -20 -10 0 10 20 30 40 IVV (kA) Net eddy current in the VV ph2 (f) Net eddy current in the VV of phase 2. Figure 6.4. Coiset currents given as inputs and RZIp outputs for both phases, plasma current and the net eddy current induced in the vacuum vessel. The desired plasma current is sustained for 20ms in phase 1 and for 100ms in phase 2 (flat-top time). 47 6.2. BREAK-DOWN RESULTS CHAPTER 6. RESULTS -40 -20 0 20 40 60 80 time (ms) 0 5 10 15 20 25 30 35 40 Ip (kA) Plasma current ph1 (a) Plasma current of phase 1. 0 50 100 150 time (ms) 0 20 40 60 80 100 Ip (kA) Plasma current ph2 (b) Plasma current of phase 2, indicating the discrete points of each of the waveforms. -40 -20 0 20 40 60 80 t (ms) -1 0 1 2 3 4 5 I (kA/ms) (plasma current) ph1 (c) Plasma current gradient of phase 1. 0 50 100 150 t (ms) -2 -1 0 1 2 3 4 I (kA/ms) (plasma current) ph2 (d) Plasma current gradient of phase 2. Figure 6.5. Plasma current and plasma current gradient for both phases. phase 2 since this phase has larger Bϕ. 6.2.2 Connection length In this section, the connection length computed using both field line tracing and the empirical formula will be showed, as well as the potential gained by the electrons along their path. Figure 6.7 shows the connection length obtained by field line tracing, the electric potential gained by the electrons and the grid used for the integration. The grid contains a meshgrid of 20x20 points, and the points located where the coils are have been removed. The resolution of the grid is 0.034m in the R axis and 0.0842 in the Z axis, giving a resolution relative to the VV size of 5.26% in both axis. The first thing to note considering the connection length Lplot is that it is not symmetric with respect to the Z= 0 plane . This is a consequence of the poloidal field, its components, BRand BZ. Figure 6.8 shows the vertical and radial field a t=0ms for phase 1 (it displays the 48 6.2. BREAK-DOWN RESULTS CHAPTER 6. RESULTS (a) Toroidal field of phase 1. BT(RGeo) = 0.1T. (b) Toroidal field of phase 2. BT(RGeo) = 0.3T. (c) Poloidal field of phase 1. IVV(t= 0) = 1.6kA. (d) Poloidal field of phase 2. IVV(t= 0) = 2.5kA. (e) Lloyd’s criteria of phase 1. (f) Lloyd’s criteria of phase 2. Figure 6.6. Magnetic fields and Lloyd’s criteria at t=0ms. Due to the symmetry with respect to Z=0, only the upper portion is shown. The poloidal field null region is dashed in green. Note in this plot the small gaps on the VV due to rounding issues can be seen. 49 6.2. BREAK-DOWN RESULTS CHAPTER 6. RESULTS (a) Connection length of phase 1. (b) Connection length of phase 2. (c) U/Vloop of phase 1. (d) U/Vloop of phase 2. 0.2 0.4 0.6 0.8 1 R (m) -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 Z (m) Cross-section (e) Grid for the line tracing. It is a 20x20 grid, and the points where the coils are located in have been removed. Figure 6.7. Connection length by line tracing and electric potential for SMART at t=0ms. The computational grid is also shown. The poloidal field null region is dashed in green. Note the connection length plot displays a region near the lower PF2 coil with high values, the lines wind up around the PF2 coil without colliding with it. This region do not appears in the electric potential plot confirming it is not relevant for the break-down. 50 6.2. BREAK-DOWN RESULTS CHAPTER 6. RESULTS same behaviour for phase 2). The Sol coil produces negative vertical field in the VV, and the PF and Div coils try to compensate it by creating positive vertical field. The resultant vertical field is mostly negative in the whole VV, but its magnitude has reduced. This explain the fact that the Lplot is not symmetric, because the lines starting in the upper portion of the VV have on average higher connection length. Hence, the highest values of Lare located at the upper portion of the VV. Higher values are obtained near the inner wall since there the toroidal field has its greatest values (Bϕ∝1/R), making Bθless relevant, resulting in lines that mostly goes in the toroidal direction, with little deviation. (a) Radial field. (b) Vertical field. Figure 6.8. Radial and vertical magnetic fields at t=0ms. The vertical component of the Sol field is negative, and the PF and Div coils try to null it by creating positive vertical field. The resultant vertical is mostly negative in the whole VV, but its magnitude has been reduced. The radial field is anti-symmetric and the vertical field symmetric with respect to the Z= 0 plane. Nevertheless, the highest values appear in the surroundings of the lower PF2 coil (Z= −0.6m), because of the inhomogeneity of the field at the surroundings of a coil. Further from the coils, the field is more homogeneous, but at the surroundings on a coil, the field is mostly the field of the coil. Plots of several field lines are shown on figure 6.9, showing that some lines in the surroundings of the lower PF2 coil wind around them without colliding with them (the magnetic lines are twisted around the lower PF2 coil in a similar way the magnetic lines are twisted in a tokamak plasma). This figure also shows that, although the upper portion of the VV in general have higher connection lengths due to the fact that the vertical field is negative, if the lines start too high, they end up colliding rapidly with the upper coils or the upper VV wall. Considering the U/Vloop plot, we see that the weird region surrounding the lower PF2 coil has low value of U/Vloop in comparison to the value of the upper portion of the VV near the 51 6.3. DISCUSSION CHAPTER 6. RESULTS (10−3.6713T) in the null region, which is very similar to GlobusM value. However, although the agreement between our values and VEST and GlobuM values seems good enough, there is a source of imprecision in our simulations, and is the fact that RZIp do not run in a selfconsistent manner with the eddy currents, as commented on section 6.2.1, so the null currents do not take into account the eddy currents, and the addition of the eddy currents worsen the field, so that the lowest poloidal field is not in the field null region. Another source of error is that our simulations use the ideal current waveforms, while more realistic simulations, such as VEST simulations, also simulate the power supplys of the coils. As a conclusion for this magnetic field comparison, other effects such as errors in the coil windings or magnetic material materials surrounding the VV could also alter the magnetic field. Considering now the connection length Lby line tracing, figure 6.14 include a comparison with NSTX and phase 2, because NSTX has Bϕ= 0.3T [42]. The values of both phase 2 and NSTX are about km long, but the appearance of both plots are not very similar. This difference is caused by the fact that SMART has most of its coils inside the VV, while NSTX has its coils outside, as can be seen on (d). To prove this point, a Lplot of an older SMART configuration with all the coils outside (phase 1) is included (c), showing similarities with the NSTX plot; both plots displays high Larm-like regions, pointing upward and outward (to the outer wall). However, SMART do not have an arm pointing downward while NSTX do, and this arm seems to be caused by a lower PF coil of NSTX, PF1B, which do not appears in the upper portion of the NSTX cross-section (i.e., NSTX cross-section do not displays symmetry with respect to the Z= 0 plane). The connection length computed by the empirical formula is roughly an order of magnitude lower than the connection length obtained by field line tracing. Differences larger than a factor of two were predicted in [23] and differences about an order of magnitude, as obtained in this work, were obtained in ITER simulations [25]. Regarding Paschen’s curve and avalanche time, SMART loop voltage seem a bit low in comparison to GlobusM, but since the loop voltage is strongly dependent on the configuration of the machine, which is different for each machine, is not such a determinant factor. Instead, studying other variables such as the avalanche time would give more information, and in our case the avalanche time confirm that the break-down phase could be finished in a short enough time interval, going from 1 to 4ms for pre-fill pressures about 10−4Torr. Data reconstruction in a VEST discharge using ECRH assistance shows that the plasma has already formed and moved to the outer wall 2ms after the ramp-down of the Sol, so the avalanche time must be ≤2ms [26]. This time interval lies between the range of values obtained for SMART, with the difference that the SMART results have been obtained with an ohmic startup, without any 58 6.3. DISCUSSION CHAPTER 6. RESULTS (a) Poloidal field of SMART phase 1. (b) Poloidal field of VEST. 1G= 10−4T. (c) Poloidal field of SMART phase 2. (d) Poloidal field of GlobusM. (e) Lloyd’s criteria of SMART phase 1. (f) Lloyd’s criteria of VEST. Figure 6.13. Comparison of the poloidal field and Lloyd’s criteria with VEST [26] and GlobusM [41] (both simulated data) at the exact time when the loop voltage is induced. Phase 1 is used to compare with VEST because it the same Bϕ= 0.1T. Eϕ= 1.12V/m for VEST and 1.06V/m for phase 1. Bϕ<0.62T for GlobusM, so it is compared with phase 2, which has 0.3T. 59 6.3. DISCUSSION CHAPTER 6. RESULTS assistance. However, SMART will also have ECRH assistance for the break-down just in case it is needed to ionize the pre-fill gas. Nevertheless, the loop voltage calculated here is the loop voltage induced by only the Sol, excluding the effect of the eddy currents, which will decrease the loop voltage according to Lenz’s law. A more precise simulation should compute the loop voltage in a self-consistent manner with the eddy currents. Finally, the simulations carried out here, and in most of the bibliography regarding breakdown studies, are done ignoring the time dependence of all the variables such as the connection length or the poloidal field, so they have to be treated as a first approximation to the problem of the tokamak start-up (to deduce (3.7), the time dependence of the loss rate was neglected). Also, the plasma current needs to be considered since shortly after the loop voltage is induced, the plasma current will be high enough to create non-negligible field as compared to the structure’s field. For deeper studies of the start-up phase, specific codes have been developed to simulate all the stages of the start-up phase such as DYON [21] or DINA, which consider the time dependence of the loss rate and the role of the plasma current in the magnetic fields. 60 6.3. DISCUSSION CHAPTER 6. RESULTS (a) Connection length of SMART phase 2. (b) Connection length of NSTX. (c) Connection length of an older SMART phase 1 with the coils outside (S1-014). (d) Cross-section of NSTX. Figure 6.14. Comparison of the connection length plot with NSTX [28,43] (simulated data). Bϕ= 0.3T for NSTX and for phase 2. 61 Chapter 7 Summary and conclusions The break-down phase of the future SMART reactor of the University of Seville has been modelled numerically using the Fiesta toolbox. The fundamentals of tokamak physics and tokamak start-up have been reviewed to understand the foundation of this work. The work of this thesis summarizes months of working in the design of the SMART reactor along with the SMART team of the PSFT group, composed of physicists and engineers students and researchers, showing the most updated scenarios for the initial operational phase of SMART and a future upgrade. The coil currents have been optimized to achieve the same target equilibrium on both phases and to ensure the break-down of the pre-fill gas, hydrogen, without any assistance method. For pre-fill pressures about 10−4Torr, the gas will last between 2 and 4 ms to break-down. The eddy currents have been included in the calculations, since they play an important role in the break-down phase altering the poloidal magnetic field inside the vacuum vessel (VV). Poloidal field values around 1.1 and 2.1·10−4T are obtained in the poloidal field null region. The connection length has been computed by field line tracing and by using an empirical formula, giving the first method about an order of magnitude greater values than the second one. The electric potential gained by an electron following the field lines has also been calculated to estimate where the gas will break-down and turns into a plasma. The gas will break-down in the upper portion of the VV, near the inner VV wall. This work presents the first approach to study the break-down phase of the tokamak startup. Further studies such as the use of specific codes to model the entire tokamak start-up phase are required to ensure proper plasma initiation. Also, a deeper understanding of the Fiesta toolbox is needed for self-consistent simulations. For more realistic simulations, the power supply’s behaviour should be included in the simulations, as well as a more realistic VV. 62 Chapter 8 Bibliography [1] G.R. Satchler. Introduction to nuclear reactions. Macmillan, 2008. [2] Kenneth S. Krane. Introductory nuclear physics. John Wiley & sons, 2008. [3] K. Miyamoto. Plasma physics for controlled fusion. Springer, 2004. [4] J. Wesson. Tokamaks. Oxford University Press, Oxford, 2004. [5] Francis F. Chen. Introduction to Plasma Physics and Controlled Fusion. Springer International Publishing, 2016. [6] Roger A. Freedman Hugh D. Young. University Physics, with Modern Physics. Pearson, 2020. [7] A. Pironti and M. Walker. Fusion, tokamaks, and plasma control: an introduction and tutorial. IEEE Control Systems Magazine, 25(5):30–43, 2005. [8] A. Beghi and A. Cenedese. Advances in real-time plasma boundary reconstruction: from gaps to snakes. IEEE Control Systems Magazine, 25(5):44–64, 2005. [9] M. Ariola and A. Pironti. Plasma shape control for the jet tokamak: an optimal output regulation approach. IEEE Control Systems Magazine, 25(5):65–75, 2005. [10] G. Ambrosino and R. Albanese. Magnetic control of plasma current, position, and shape in tokamaks: a survey or modeling and control approaches. IEEE Control Systems Magazine, 25(5):76–92, 2005. [11] Linjin Zheng. Advanced tokamak stability theory. Morgan & Claypool Publishers, 2015. 63 CHAPTER 8. BIBLIOGRAPHY [12] Zhe Gao. Compact magnetic confinement fusion: Spherical torus and compact torus. Matter and Radiation at Extremes, 1(3):153 – 162, 2016. [13] Y-K.M. Peng and D.J. Strickler. Features of spherical torus plasmas. Nuclear Fusion, 26(6):769–777, jun 1986. [14] A. Mancini, J. Ayllon-Guerola, S. J. Doyle, M. Agredano-Torres, D. Lopez-Aires, J. Toledo- Garrido, C. Soria-Hoyo, E. Viezzer, M. Garcia-Munoz, PF. Buxton, K.-J. Chung, J. Garcia-Dominguez, MP. Gryaznevich, J. Hidalgo-Salaverri, Y.-S. Hwang, and J. Segado- Fern´andez. Mechanical and electromagnetic design of the vacuum vessel of the smart tokamak. In 31st Symposium on Fusion Technology (SOFT-31), 2020. [15] Jefrrey P. Freidberg. Ideal magnetohydrodynamics. Plenum Press, 1987. [16] T C Luce. An analytic functional form for characterization and generation of axisymmetric plasma boundaries. Plasma Physics and Controlled Fusion, 55(9):095009, jul 2013. [17] Y.-K. M. Peng. The physics of spherical torus plasmas. Physics of Plasmas, 7(5):1681– 1692, 2000. [18] D. Mueller. The physics of tokamak start-up. Physics of Plasmas, 20(5):058101, 2013. [19] P.C. de Vries and Y. Gribov. Iter breakdown and plasma initiation revisited. Nuclear Fusion, 59(9):096043, aug 2019. [20] Joyeeta Sinha. Plasma breakdwon and current formation in single core and doublet configurations on TCV. PhD thesis, ´ Ecole Polytechnique F´ed´erale de Laussane, 2017. [21] Hyun-Tae Kim. Physics and computational simulatios of plasma burn-through for tokamak start-up. PhD thesis, Imperial college of Science Technology and Medicine, London, 2013. [22] B Lloyd, P G Carolan, and C D Warrick. Ecrh-assisted start-up in ITER. Plasma Physics and Controlled Fusion, 38(9):1627–1643, sep 1996. [23] B. Lloyd, G.L. Jackson, T.S. Taylor, E.A. Lazarus, T.C. Luce, and R. Prater. Low voltage ohmic and electron cyclotron heating assisted startup in diii-d. Nuclear Fusion, 31(11):2031–2053, nov 1991. [24] ITER Physics Expert Group on Disruptions, Plasma Control, ITER Physics Expert Group on Energetic Particles, Heating MHD, Current Drive, ITER Physics Expert Group on Di- 64 CHAPTER 8. BIBLIOGRAPHY agnostics, and ITER Physics Basis Editors. Chapter 8: Plasma operation and control. Nuclear Fusion, 39(12):2577–2625, dec 1999. [25] Y Gribov, D Humphreys, K Kajiwara, E.A Lazarus, J.B Lister, T Ozeki, A Portone, M Shimada, A.C.C Sips, and J.C Wesley. Chapter 8: Plasma operation and control. Nuclear Fusion, 47(6):S385–S403, jun 2007. [26] YoungHwa An, Jeongwon Lee, HyunYeong Lee, JongGab Jo, Bong-Ki Jung, Kyoung-Jae Chung, Young-Gi Kim, Jungmin Jo, Jeong hun Yang, Yong-Su Na, T.S. Hahm, and Y.S. Hwang. Plasma start-up design and first plasma experiment in vest. Fusion Engineering and Design, 96-97:274 – 280, 2015. Proceedings of the 28th Symposium On Fusion Technology (SOFT-28). [27] E.A Lazarus, A.W Hyatt, G.L Jackson, and D.A Humphreys. Using a multipole expansion for startup in the diii-d tokamak. Nuclear Fusion, 38(7):1083–1096, jul 1998. [28] K.C. Hammond, R. Raman, and F.A. Volpe. Application of townsend avalanche theory to tokamak startup by coaxial helicity injection. Nuclear Fusion, 58(1):016013, nov 2017. [29] David J. Griffiths. Introduction to electrodynamics. Prentice Hall, New Jersey, 1999. [30] L. Pangione, G. McArdle, and J. Storrs. New magnetic real time shape control for mast. Fusion Engineering and Design, 88(6):1087 – 1090, 2013. Proceedings of the 27th Symposium On Fusion Technology (SOFT-27); Li`ege, Belgium, September 24-28, 2012. [31] M J Windridge, G Cunningham, T C Hender, R Khayrutdinov, and V Lukash. Nonlinear instability at large vertical displacements in the mast tokamak. Plasma Physics and Controlled Fusion, 53(3):035018, feb 2011. [32] L.L. Lao, H. St. John, R.D. Stambaugh, A.G. Kellman, and W. Pfeiffer. Reconstruction of current profile parameters and plasma shapes in tokamaks. Nuclear Fusion, 25(11):1611– 1622, nov 1985. [33] Alfredo Pironti Marco Ariola. Magnetic control of tokamak plasmas. Springer, 2008. [34] Atul Stefan Sharma. Tokamak modelling and control. PhD thesis, Imperial college of Science Technology and Medicine, London, 2002. [35] J.B. Lister, A. Sharma, D.J.N. Limebeer, Y. Nakamura, J.P. Wainwright, and R. Yoshino. Plasma equilibrium response modelling and validation on jt-60u. Nuclear Fusion, 42(6):708–724, jun 2002. 65 CHAPTER 8. BIBLIOGRAPHY [36] S. J. Doyle, A. Mancini, M. Agredano-Torres, J. L. Garcia-Sanchez, J. M. Ayllon-Guerola, D. Lopez-Aires, M. Garcia-Munoz, E. Viezzer, C. Soria-Hoyo, P. F. Buxton, M. P. Gryaznevich, Y.-S. Hwang, and K.-J. Chung. Equilibrium design for the smart tokamak. In 31st Symposium on Fusion Technology (SOFT-31), 2020. [37] Peter Buxton. Merging/compression start-up in st40: Analysis of first experimental results. In 45th EPS Conference on Plasma Physics, 07 2018. [38] K.J. Chung, Y.H. An, B.K. Jung, H.Y. Lee, J.J. Dang, J.W. Lee, J. Yang, J.G. Jo, D.H. Choi, Y.G. Kim, Y.S. Na, and Y.S. Hwang. Initial plasma start-up using partial solenoid coils in versatile experiment spherical torus (vest). Fusion Engineering and Design, 88(6):787 – 790, 2013. Proceedings of the 27th Symposium On Fusion Technology (SOFT- 27); Li`ege, Belgium, September 24-28, 2012. [39] V Mertens, Guenter Haas, V. Rohde, and ASDEX Team. Hydrogen gas balance in asdex upgrade with div iib. In Koch, R.; Lebedev, S.: 30th EPS Conference on Plasma Physics and Controlled Fusion, European Physical Society (2003), 01 2003. [40] X. Song, X.R. Duan, X.M. Song, G.Y. Zheng, S. Wang, B. Li, X.Y. Bai, Sh.D. Song, C. Wang, and J. Sun. Experimental results of plasma breakdown and flux optimization on hl-2a tokamak. Fusion Engineering and Design, 125:195 – 198, 2017. [41] V.K Gusev, T.A Burtseva, A.V Dech, G.A Gavrilov, V.E Golant, S.V Krikunov, R.G Levin, V.B Minaev, A.B Mineev, O.A Minyaev, E.E Mukhin, A.N Novokhatskii, Yu.V Petrov, E.N Rumyantsev, N.V Sakharov, and V.M Sharapov. Plasma formation and first oh experiments in the globus-m tokamak. Nuclear Fusion, 41(7):919–925, jul 2001. [42] Stanley M. Kaye, Masayuki Ono, Yueng-Kay Martin Peng, Donald B. Batchelor, Mark D. Carter, Wonho Choe, Robert Goldston, Yong-Seok Hwang, E. Fred Jaeger, Thomas R. Jarboe, Stephen Jardin, David Johnson, Robert Kaita, Charles Kessel, Henry Kugel, Rajesh Maingi, Richard Majeski, Janhardan Manickam, Jonathan Menard, David R. Mikkelsen, David J. Orvis, Brian A. Nelson, Franco Paoletti, Neil Pomphrey, Gregory Rewoldt, Steven Sabbagh, Dennis J. Strickler, Edmund Synakowski, and James R. Wilson. Physics design of the national spherical torus experiment. Fusion Technology, 36(1):16–37, 1999. [43] D. A. Gates, J. E. Menard, and R. J. Marsala. Vessel eddy current measurement for the national spherical torus experiment. Review of Scientific Instruments, 75(12):5090–5093, 2004. 66 CHAPTER 8. BIBLIOGRAPHY [44] Norman S. Nise. Control system engineering. Wiley, 2019. 67