Evolution of Heavy Ion Beam Probing from the Origins to Study of Symmetric Structures in Fusion Plasmas

The overview discusses development of the unique fusion plasma diagnostics—Heavy Ion Beam Probing (HIBP) in application to toroidal magnetic plasma devices. The basis of the HIBP measurements of the plasma electric potential and processing of experimental data are considered. Diagnostic systems for probing plasma in tokamaks TM-4, TJ-1, TUMAN-3M and T-10, stellarators WEGA, TJ-II and Uragan-2M are presented. Promising results of the HIBP projects for various existing modern machines, such as TCV, TCABR, MAST, COMPASS, GLOBUS-M2, T-15 MD and W7-X and the international fusion tokamak reactor ITER are given. Results from two machines with similar size and plasma parameters, but with different types of the magnetic con-figuration: axisymmetric tokamak T-10 and helically symmetric stellarator TJ-II are compared. The results of studies of stationary potential profiles and oscillations in the form of quasimonochromatic and broadband fluctuations, turbulent particle flux, fluctuations of density and poloidal magnetic field are presented. The properties of symmetric structures—zonal flows and geodesic acoustic modes of plasma oscillations as well as Alfvén Eigenmodes excited by fast particles from neutral beam injection heating are described. General trends in the behavior of electric potential and turbulence in magnetized fusion plasmas are revealed.


Preface
In recent decades, diagnostics of high-temperature plasma has become an independent scientific discipline that has grown out of the physics of fusion plasma. Fusion plasma is a very unusual area of science. Unlike other areas of physics that study the properties of objects existing on Earth, plasma physics first creates an object of its research in terrestrial conditions-plasma with thermonuclear fusion temperatures higher than the temperature in the Sun core. Superconducting magnets, creating magnetic fields for plasma confinement, operate near the absolute minimum temperature, in conditions of space cold. The vacuum conditions in fusion machines for magnetic confinement of plasma are approaching to the conditions of cosmic space. The creation of such a complex and un-stable object as a fusion plasma in a terrestrial laboratory is in itself a formidable physical and technical issue. An equally difficult issue is plasma diagnostics-the development of tools and methods for studying plasma, reliable techniques and tools for measuring its properties and parameters, such as the temperature of its electrons and ions, density, stored energy, effective charge, and many others. The listed plasma parameters are not measured directly, which is obviously impossible in such a hot medium, but by their indirect signs, such as, for example, plasma radiation in different energy ranges. From the point of view of mathematics, the determination of the required plasma parameters from the quantities available for measurement is a classic example of inverse problems that require the development of special regularization methods for their solution [1,2]. In the last decades of the 20th

Introduction
The study of electric fields is one of the urgent problems of modern physics of fusion plasma. These studies are related to both the fundamental problem-the emergence of an electric field in a quasineutral plasma [3], and an applied problem-the understanding of the electric field effects on the plasma transport processes in closed magnetic configurations [4]. It is known that the anomalous mechanisms linked with plasma turbulence are dominated in the particle and energy transport across confining electric field [5]. Currently, the hypothesis that the plasma turbulence and anomalous transport are stabilized by shear (radial inhomogeneity) of poloidal plasma rotation in crossed radial electric E r and toroidal magnetic B t fields is generally accepted in magnetic confinement studies [6]. However, this hypothesis has both agreements and disagreements with experiment, so, the effect of E r to plasma confinement still remains unclear.
At present, it is believed that the problems associated with the electric field in tokamaks are becoming the principal ones in fusion research, and that for the success of the ITER project, it is necessary to change from studying the role of the electric field to controlling the anomalous transport with electric fields [7]. Measurement of the electric field in the plasma of modern fusion machines is a complex experimental issue [8]. At the plasma edge, the electric potential is measured with Langmuir probes. In the hot plasma core, the electric field, as a rule, is not measured directly, but it is determined from the measured plasma rotation velocity using spectroscopy or correlation reflectometry. The only direct method to find the electric potential in the plasma interior is the Heavy Ion Beam Probe (HIBP) diagnostics [9].
Diagnostic techniques for toroidal plasmas always have limitations concerning the access to the plasma. These constrains are dictated by the complexity of the toroidal devices and by the trend to reduce the ripple of the confining magnetic field by an increase in the number of toroidal field coils. This makes the diagnostic ports narrow in toroidal direction. Finally, most of the diagnostics observe just a part of the plasma vertical cross-section. They have to assume the symmetry of the measured parameter to complete the whole picture of its spatial distribution. The typical example of such an approach is the multi-chord plasma density measurement with a microwave interferometer as a function of the plasma radius that then assumes a poloidal symmetry of the density distribution.
On top of that, there are theoretical predictions of various asymmetries for plasma parameters. For example, there could be the so-called ballooning effect, which deforms the contour lines of plasma parameters at the Low-Field Side (LFS) of the torus vertical cross-section, breaking the in-out symmetry. Gyrokinetic simulations predict a strong increase of the broadband electrostatic turbulence in the LFS in respect to the High Field Side (HFS). Therefore, the issue of symmetry is an open one in many aspects of plasma studies. Resolving this issue is highly demanding for the diagnostic techniques.
This overview describes the HIBP diagnostics of toroidal plasma and its development over the past four decades for studying various phenomena in plasma including issues of their symmetry.

Heavy Ion Beam Probe-A Tool for Measuring Electrical Potential and Plasma Turbulence
The Heavy Ion Beam Probe (HIBP) is a unique method for plasma studies in fusion machines [10]. HIBP was first implemented in the late 1960s by Hickok and Jobs on the arc discharge and later on the ST tokamak [11]. For the first time in Russia (USSR), HIBP was implemented in the early 1980s by Krupnik and Nedzelskiy at the TM-4 tokamak [12,13]. Currently, HIBP is the only non-disturbing method for direct measurements of the electric potential in a hot plasma core [14]. In this paper, we consider the mathematical aspects of diagnostics, equipment and methods for measuring plasma parameters available for the modern HIBP [15].

Physical Principles of HIBP Measurements
To probe the plasma with the beam of heavy ions, the beam is injected across the confining magnetic field. When the beam particles fly through the plasma, some of them collide with plasma particles (mainly electrons) and lose one or more electrons. As a result, a fan of secondary ionized particles with higher charges is formed. The HIBP diagnostic scheme and an example of its realization in the T-10 tokamak are shown in Figure 1. This overview describes the HIBP diagnostics of toroidal plasma and its development over the past four decades for studying various phenomena in plasma including issues of their symmetry.

Heavy Ion Beam Probe-A Tool for Measuring Electrical Potential and Plasma Turbulence
The Heavy Ion Beam Probe (HIBP) is a unique method for plasma studies in fusion machines [10]. HIBP was first implemented in the late 1960s by Hickok and Jobs on the arc discharge and later on the ST tokamak [11]. For the first time in Russia (USSR), HIBP was implemented in the early 1980s by Krupnik and Nedzelskiy at the TM-4 tokamak [12,13]. Currently, HIBP is the only non-disturbing method for direct measurements of the electric potential in a hot plasma core [14]. In this paper, we consider the mathematical aspects of diagnostics, equipment and methods for measuring plasma parameters available for the modern HIBP [15].

Physical Principles of HIBP measurements
To probe the plasma with the beam of heavy ions, the beam is injected across the confining magnetic field. When the beam particles fly through the plasma, some of them collide with plasma particles (mainly electrons) and lose one or more electrons. As a result, a fan of secondary ionized particles with higher charges is formed. The HIBP diagnostic scheme and an example of its realization in the T-10 tokamak are shown in Figure 1. A small-aperture detector is usually placed outside the magnetic field that allows the secondary ion beam (part of the fan), formed at the secondary ionization point (sample volume or SV) on the probe beam trajectory, to enter the detector. The detected secondary A small-aperture detector is usually placed outside the magnetic field that allows the secondary ion beam (part of the fan), formed at the secondary ionization point (sample volume or SV) on the probe beam trajectory, to enter the detector. The detected secondary ions carry information about the plasma parameters in the SV. The spatial resolution of the method λ is determined by the dimensions of the sample volume and mainly depends on the size of the aperture and the position of the detector. In real experiments on different machines λ = 0.2-2 cm. HIBP allows us to simultaneously measure several plasma parameters: electric potential ϕ, electron density n e , and magnetic potential A (or poloidal magnetic field B p or current density j). The position of the sample volume can be moved along the plasma cross section by varying the energy of the probing beam E b or the injection angle of the beam into the plasma α. Thus, the set of points observed by the detector forms a two-coordinate (E b , α) grid called the detector grid.
The trajectories of primary and secondary ions in a low-current tokamak lie in the meridional plane between the confining field coils or near it. For tokamaks with a high current density and for stellarators, the trajectories of the probing particles are shifted in the toroidal direction, and the detector grid becomes three-dimensional. This essential feature of the motion of probing particles should be taken into account during development projects of placing diagnostics on new machines installations, and in the analysis of measurement results.
Thus, a specific feature of HIBP diagnostics is that the values of plasma parameters are determined from measuring the characteristics of the secondary beam: intensity, energy, etc., while the position of the measurement point can only be found from calculations of the trajectories of probing particles. Therefore, before placing HIBP diagnostics on a new machine installation, not only the development of diagnostic equipment is required, but also preliminary calculations.
The position of the sample volume in the meridional plane depends on a large number of probing parameters, including: (i) geometric parameters: the coordinates of the injection x i , y i and detection points x D , y D , and the injection angle α (ii) the physical parameters: the toroidal magnetic field B t , the energy E b , the mass m, and the electric charge q of the probing particles.
Probing particles move in the magnetic field of the machine along a Larmor circle with a radius: It is necessary to determine not the value of R L itself, but the values of the parameters on the right hand side of Equation (1), which are subject to independent restrictions. The equality R L (q, m, E b , B t ) = const reduces the number of physical parameters from 4 to 3. So, for 2D optimization in the meridional plane, there are eight parameters: x I , y I , x D , y D , α, B t , E b , q. The need to pass the trajectories of particles through the ports of the vacuum chamber and the structural elements of the machine impose serious constraints on these parameters. These constraints are related to each other because the allowable tolerances for each parameter depend on the others. Such relations defy analytical description due to the complexity of the trajectory behavior. The configuration of the vacuum chamber in the meridional plane approximately determines their tolerances. The problem of optimizing experimental conditions includes the choice of optimal probing parameter values in the sense of the following goal functions:

1.
Pass the beam through the existing vacuum vessel ports; 2.
Find the detector line from the core to the edge of the plasma; 3.
Find a detector grid covering the maximum part of the plasma cross section; 4.
Optimize the range of beam energies.
In high-current tokamak discharges, as well as in stellarators, the magnetic field has a noticeable poloidal component B p , which shifts the orbits of the probing particles out of the meridional plane and transforms them into spatial curves [26]. In this case, additional parameters appear that affect the sample volume. These are the toroidal coordinates of the injector and detector, Z I and Z D , and the injection angle between the initial velocity vector and the meridional plane β. This gives us a 12-dimensional optimization problem with implicit connections between the limits of parameters and with non-formalized goal functions. Figure 2 shows an example of solving this problem by 'the shooting' method in T-10. Examples of solving the optimization problem in other machines are given in Sections 7 and 8. The first goal function was achieved in all considered machines, as the chosen probing schemes ensured the passage of trajectories through the vacuum vessel ports and allowed measurements. Maximization of the second and third goal functions provides the maximum possible radial and angular sizes of the observation area in plasma. The fourth function is not always so essential, because the energy range of modern accelerators is quite large. However, for large machines, such as the W-7X stellarator and the constructed ITER tokamak reactor, this task is of great importance. As will be shown in Section 7, the chosen injection schemes allow one to connect the core and the edge of the plasma with a detector line using the acceptable beam energy.
vector and the meridional plane β. This gives us a 12-dimensional optimization problem with implicit connections between the limits of parameters and with non-formalized goal functions. Figure 2 shows an example of solving this problem by 'the shooting' method in T-10. Examples of solving the optimization problem in other machines are given in Sections 7 and 8. The first goal function was achieved in all considered machines, as the chosen probing schemes ensured the passage of trajectories through the vacuum vessel ports and allowed measurements. Maximization of the second and third goal functions provides the maximum possible radial and angular sizes of the observation area in plasma. The fourth function is not always so essential, because the energy range of modern accelerators is quite large. However, for large machines, such as the W-7X stellarator and the constructed ITER tokamak reactor, this task is of great importance. As will be shown in Section 7, the chosen injection schemes allow one to connect the core and the edge of the plasma with a detector line using the acceptable beam energy.

Determination the Spatial Distribution of Electric Potential
The principle of HIBP measurements is based on the potentiality of the electric field, the conservation of the total energy of the charged particles in the plasma electric field ( Figure 3). The probing beam, hereinafter referred to as the primary beam, enters the plasma with an initial energy Eb. At the ionization point, which is the plasma region of interest (SV), a particle of the primary beam loses an electron with potential energy which is acquired by the secondary ionized probe ion. The total energy of secondary ions leaving the plasma is . Therefore, the local potential in SV is equal to the energy difference:

Determination the Spatial Distribution of Electric Potential
The principle of HIBP measurements is based on the potentiality of the electric field, the conservation of the total energy of the charged particles in the plasma electric field ( Figure 3). The probing beam, hereinafter referred to as the primary beam, enters the plasma with an initial energy E b . At the ionization point, which is the plasma region of interest (SV), a particle of the primary beam loses an electron with potential energy −eϕ SV pl which is acquired by the secondary ionized probe ion. The total energy of secondary ions leaving the plasma is E d = E b + eϕ SV pl . Therefore, the local potential in SV is equal to the energy difference: Symmetry 2021, 13, x FOR PEER REVIEW 6 of 49 Thus, the problem of measuring the plasma electric potential does not pose serious mathematical difficulties. However, from the point of view of conducting physical measurements it is a very difficult task, as the dimensions of modern plasma machines (scale of 15 m) and the magnitude of the magnetic field (scale of 1-3 T) impose high requirements on the Larmor radius of probing particles, which must exceed the plasma size. According to (1), this requires the use of probing ions with the minimum charge number (q Thus, the problem of measuring the plasma electric potential does not pose serious mathematical difficulties. However, from the point of view of conducting physical Symmetry 2021, 13, 1367 6 of 46 measurements it is a very difficult task, as the dimensions of modern plasma machines (scale of 1-5 m) and the magnitude of the magnetic field (scale of 1-3 T) impose high requirements on the Larmor radius of probing particles, which must exceed the plasma size. According to (1), this requires the use of probing ions with the minimum charge number (q = 1), maximum possible mass (Cs + , Tl + , Au + ) and energies (on a scale of several hundred keV-several MeV). On the other hand, the measured plasma potential has the scale of the electron temperature. In the plasma core, it ranges from several hundred eV to keV. At the plasma edge, it has a scale from several eVs to hundreds of eV. Thus, according to (2), the measurement of the potential claims the computation of a small difference of two large quantities. Such measurements require high accuracy of each quantity or highly stabilized high-voltage sources and the use of precision instruments to measure the energy of the secondary beam.

Determination of the Spatial Distribution of Plasma Density
The local values of the plasma density n(l) can be found from the relation: where i 1 and i 2 are the intensities of the primary and secondary beams in the SV; σ = συ e /υ b is the effective ionization cross section by electron impact, <συ e > is the ionization rate averaged over the Maxwellian distribution of the electron velocity υ e , itself a function of the electron temperature T e ; υ b is the initial velocity of the probing beam, υ e > υ b ; l = l(x, y) is the coordinate on the detector line; and λ(l(E)) is the longitudinal dimension of SV, the length of the arc along the primary trajectory, from which the secondary ions reach the detector. The electron density n(l) is averaged over the length λ.
If the density is low and the dimensions of the machine are small, then we can assume that the intensities i 1 and i 2 are the intensities of the primary beam at the exit from the ion injector and the secondary beam at the detector. In this case, Relation (3) can be directly used to determine the local values of the plasma density. However, in modern fusion machines installations, both the primary and secondary beams are significantly attenuated, when passing through the plasma. Methods for determining plasma density using HIBP in the case of strong attenuation are described in [27,28], and they were used for interpretation of experimental data from the TM-4 tokamak.
The expression for the current of secondary ions to the detector I t , taking into account the attenuation of along the beam trajectory, caused by elementary processes, including collisions of beam particles with plasma particles, can be written as follows: where I b is the ion beam current leaving the injector, n is the electron density, σ 12 and σ 23 are the known effective cross sections for collisional ionization of the probing beam (e.g., Cs in the cases of TM-4 and TJ-II) by plasma electrons.
The unknown function in this nonlinear integro-functional equation is the plasma electron density n(l); functional dependence I t (l) is determined experimentally; I b is a known constant. The integration lines Il and lD, and the quantities λ(l) are found by calculating the trajectory and the detector line.
Equation (4) is an example of a nonlinear inverse problem with a nonclassical operator structure: the integral of the desired function is in the exponent. The limits of integration are variable; they depend on the coordinates of a point on the detector line l = l (x, y). This equation may be solved by a quasi-solution method: parametrization the desired function and finding the optimal values of the parameters by minimizing the discrepancy functional. The application of this method to the analysis of experimental data from TM-4 is shown in Figure 4. We see that the reconstructed density profile is consistent with the data of other diagnostics, so the method works well for TM-4 conditions. More recently, the same technique was applied for the TJ-II data. The results show that the reconstructed HIBP density profile lies within the error bars of the Thomson scattering data [29]. electron density n(l); functional dependence It(l) is determined experimentally; Ib is a known constant. The integration lines Il and lD, and the quantities λ(l) are found by calculating the trajectory and the detector line.
Equation (4) is an example of a nonlinear inverse problem with a nonclassical operator structure: the integral of the desired function is in the exponent. The limits of integration are variable; they depend on the coordinates of a point on the detector line l = l (x, y). This equation may be solved by a quasi-solution method: parametrization the desired function and finding the optimal values of the parameters by minimizing the discrepancy functional. The application of this method to the analysis of experimental data from TM-4 is shown in Figure 4. We see that the reconstructed density profile is consistent with the data of other diagnostics, so the method works well for TM-4 conditions. More recently, the same technique was applied for the TJ-II data. The results show that the reconstructed HIBP density profile lies within the error bars of the Thomson scattering data [29].

Determination of the Plasma Magnetic Potential (Field of Plasma Current)
The plasma current profile mainly determines the equilibrium of the plasma column, its stability, energy balance, and the characteristics of fusion plasma. At present, for measurements of the poloidal magnetic field, there is a very limited set of diagnostics [30], which do not have a high temporal and spatial resolution. Therefore, the problem of creating diagnostics for poloidal fields in toroidal plasma is still urgent. In this Subsection, we consider a method for determining the plasma current field from measurements of the shift of the probing beam caused by this field [31,32].
Let us consider the simplest model of a symmetric torus (without ripple), in which the toroidal magnetic moment is conserved for a probing charged particle. After integration the Lagrangian for the probing particle along the trajectories, we obtain the following expression for the desired function A, the toroidal component of the magnetic potential:

Determination of the Plasma Magnetic Potential (Field of Plasma Current)
The plasma current profile mainly determines the equilibrium of the plasma column, its stability, energy balance, and the characteristics of fusion plasma. At present, for measurements of the poloidal magnetic field, there is a very limited set of diagnostics [30], which do not have a high temporal and spatial resolution. Therefore, the problem of creating diagnostics for poloidal fields in toroidal plasma is still urgent. In this Subsection, we consider a method for determining the plasma current field from measurements of the shift of the probing beam caused by this field [31,32].
Let us consider the simplest model of a symmetric torus (without ripple), in which the toroidal magnetic moment is conserved for a probing charged particle. After integration the Lagrangian for the probing particle along the trajectories, we obtain the following expression for the desired function A ϕ , the toroidal component of the magnetic potential: Here s 0 is the injection point, s i is the ionization point, and s d is the detection point, q and q + k are dimensionless charge numbers of primary and secondary beam particles. Equation (6) connects the unknown function A ϕ with the angular displacement of the beam in the detector respect to the one of injector ϕ d − ϕ 0 , that is, with the measured quantity. Equation (6) is a nonlinear Volterra-type integral equation of the second kind. As the integration limits are variable, they depend on the coordinates of the ionization point in the plasma (points on the detector line). In contrast with the equation of the first kind, solving such equations is a well-posed problem. The function A ϕ is unknown only inside the plasma column, for r < a. Outside the plasma, where r > a, A ϕ ext is determined by its values inside the plasma as follows: where the integration is carried out over the entire region of current flow, and R is the distance from the volume element dV to the position, where A ϕ ext is computed. Note that the integrals in (6-7) are curvilinear integrals over 3D trajectories. The arc element ds along the particle trajectory depends on the unknown function A ϕ : ds = ds(A ϕ ), therefore, Equation (7) is a nonlinear integral equation.
The solving of this equation is based on the iteration method. Figure 5 shows the results of the numerical solution of Equation (7) for the conditions of the TUMAN-3M tokamak.
and q + k are dimensionless charge numbers of primary and secondary beam particles. Equation (6) connects the unknown function Aφ with the angular displacement of the beam in the detector respect to the one of injector φd − φ0, that is, with the measured quantity. Equation (6) is a nonlinear Volterra-type integral equation of the second kind. As the integration limits are variable, they depend on the coordinates of the ionization point in the plasma (points on the detector line). In contrast with the equation of the first kind, solving such equations is a well-posed problem. The function Aφ is unknown only inside the plasma column, for r < a. Outside the plasma, where r > a, Aφ ext is determined by its values inside the plasma as follows: where the integration is carried out over the entire region of current flow, and R is the distance from the volume element dV to the position, where Aφ ext is computed. Note that the integrals in (6-7) are curvilinear integrals over 3D trajectories. The arc element ds along the particle trajectory depends on the unknown function Aφ: ds = ds(Aφ), therefore, Equation (7) is a nonlinear integral equation. The solving of this equation is based on the iteration method. Figure 5 shows the results of the numerical solution of Equation (7) for the conditions of the TUMAN-3M tokamak.

Application of the HIBP to Fluctuation Measurements
The HIBP is a multi-purpose diagnostic tool capable to measure the local potential φ using the secondary beam energy, plasma density ne from the beam current It, and the poloidal magnetic field Bp from the toroidal angular displacement of beam φd or the linear beam shift ζ = φd R, where R is the distance from the detection point to the origin of the cylindrical coordinate system [33].

Application of the HIBP to Fluctuation Measurements
The HIBP is a multi-purpose diagnostic tool capable to measure the local potential ϕ using the secondary beam energy, plasma density n e from the beam current I t , and the poloidal magnetic field B p from the toroidal angular displacement of beam ϕ d or the linear beam shift ζ = ϕ d R, where R is the distance from the detection point to the origin of the cylindrical coordinate system [33].
High temporal resolution (up to 0.5 µs) allows measurements of fluctuations of the parameters under study. Expressions for the links between the relative fluctuations of the measured quantities and the parameters under study are following: For the variable component of the potential: where δi is the normalized difference between the beam currents across the separated detector plates, U an is the voltage applied to the analyzer, and F is its dynamic coefficient. Thus, there is a simple linear relationship between the desired value δϕ SV pl and the measured value δi(t).
If the spatial correlation length of density fluctuations is much less than the length of the beam trajectory, then the nonlocal terms in Equation (4) may be omitted due to integration of the oscillatory component δn e over the beam trajectory. In this case: For large-scale density fluctuations δn e , for which the radial correlation length is large and comparable to the plasma radius, the local data on oscillations in the sample volume can be contaminated by attenuation effects accumulating along the entire trajectory. For the experimental data, considered in this paper the large-scale modes are not typical, so the density oscillations are measured locally [34].
The oscillating component of the toroidal shift of the beam in the detector ζ d , is described by the following expression [35]: Here, ψ = A ϕ /R, P ζ inj is the toroidal magnetic moment at the injection point, known constant. Note that the integrals in (10) are taken along the total paths L 1 and L 2 from the injector to the position of the sample volume (SV), and from SV to the detector. Equation (10) shows that the local term δψ SV input the main contribution to the oscillations of beam shift. The contribution of integral terms is much smaller for oscillations with shorter correlation length (or radial wavelength) due to the path integration. It was shown in Reference [36] that the integral terms do not strongly affect the local terms from MHD oscillations, therefore HIBP provides practically local measurements of magnetic fluctuations. However, a nonlocal contribution to ζ d from oscillations with a large correlation length, comparable to the length of the beam trajectory in plasma, is also possible.
In conclusion, we note that, despite the path integral effect, HIBP provides fairly local measurements of density and magnetic field fluctuations, taking into account the conditions specified above, as well as purely local measurements of electric potential fluctuations. All local measurements are averaged over the sample volume (SV), which size determines the spatial resolution of the HIBP diagnostic. The sample volume typically looks like the inclined elliptical disk with diameter equal to the diameter of the beam, and its thickness corresponds to the entrance slit of the analyzer. The typical radial dimension of SV ranges from 0.1 to 1 cm.

Spectral Fourier Analysis of Oscillations
To analyze the oscillating component of an arbitrary signal x(t), we introduce the Fourier transform as S x (f ) and S x *(f ) as its complex conjugate. Thus the power spectral density (PSD) of oscillations is: The cross-spectral oscillation power density (CSPD) for signals x(t) and y(t) is the following product of their Fourier transforms: An example of calculating the PSD of oscillations of plasma parameters determined in a certain time domain is shown in Figure 6a. As a rule, when digitizing an experimental signal with a frequency of 1 µs, the PSD is calculated over a period of 0.25-1 ms, depending on the characteristic time of the change in the frequency of the processes under study.
following product of their Fourier transforms: An example of calculating the PSD of oscillations of plasma parameters determined in a certain time domain is shown in Figure 6a. As a rule, when digitizing an experimental signal with a frequency of 1 μs, the PSD is calculated over a period of 0.25-1 ms, depending on the characteristic time of the change in the frequency of the processes under study. The power spectra of characteristic peaks of monochromatic fluctuations corresponding to the Geodesic Acoustic Mode (GAM) [37,38] and magnetohydrodynamic (MHD) tearing mode with a poloidal mode number m = 2 are seen. In contrast to the solitary peak of the m = 2 mode, the frequency GAM peak is accompanied by a high-frequency satellite peak. Figure 7 presents the physical picture of the zonal flows as symmetric torsional plasma oscillations [37]. While zonal flow by itself, also named residual zonal flow, typically has low-to-near zero frequency, GAM as a higher frequency branch of zonal flows has a typical frequency of around 20 kHz. Both GAM and tearing modes present examples of symmetric structure in the plasma toroidal magnetic configuration. GAM is poloidally and toroidally symmetric torsional plasma oscillation with n = m = 0, where n is the toroidal mode number [39]. MHD tearing mode is a chain of magnetic islands, that have helically symmetric structure, located on magnetic surface with rational value of safety factor q = m/n in the toroidal magnetic configuration [40]. The power spectra of characteristic peaks of monochromatic fluctuations corresponding to the Geodesic Acoustic Mode (GAM) [37,38] and magnetohydrodynamic (MHD) tearing mode with a poloidal mode number m = 2 are seen. In contrast to the solitary peak of the m = 2 mode, the frequency GAM peak is accompanied by a high-frequency satellite peak. Figure 7 presents the physical picture of the zonal flows as symmetric torsional plasma oscillations [37]. While zonal flow by itself, also named residual zonal flow, typically has low-to-near zero frequency, GAM as a higher frequency branch of zonal flows has a typical frequency of around 20 kHz. Both GAM and tearing modes present examples of symmetric structure in the plasma toroidal magnetic configuration. GAM is poloidally and toroidally symmetric torsional plasma oscillation with n = m = 0, where n is the toroidal mode number [39]. MHD tearing mode is a chain of magnetic islands, that have helically symmetric structure, located on magnetic surface with rational value of safety factor q = m/n in the toroidal magnetic configuration [40]. This spectrogram, describing a process lasting 500 ms, consists of a sequence of elementary spectra, each of which is calculated in 1 ms. The potential spectrogram is dominated by the GAM peak, accompanied by a satellite peak [41]. Application of additional ECR heating in the time interval from 600 to 800 ms leads to an increase in the frequency of both the GAM and the satellite due to growth of the electron temperature [42]. This spectrogram, describing a process lasting 500 ms, consists of a sequence of elementary spectra, each of which is calculated in 1 ms. The potential spectrogram is dominated by the GAM peak, accompanied by a satellite peak [41]. Application of additional ECR heating in the time interval from 600 to 800 ms leads to an increase in the frequency of both the GAM and the satellite due to growth of the electron temperature [42].

Fourier Analysis of Coherency and Cross-Phase of Oscillations
Time evolution of the correlations between signals x(t) and y(t) are determined by spectrograms of their coherency Coh xy and cross-phase θ xy , as following: Coherence and cross-phase are modulus and argument of the complex function CPSD. Figure 8 presents an example of calculating the coherence spectrogram. Values of statistically significant coherence indicate the Alfvén eigenmodes (AEs), which is another type of MHD mode. AEs are electromagnetic oscillations, propagating along the magnetic field line with Alfvén speed and possessing helically symmetric structure with finite m and n. In the vicinity of the magnetic flux surface with rational q value, equals to m/n for a specific mode, the dispersion relation for AE contains also a GAM component, causing GAM frequency to be a lower limit for AE frequency [43,44]. It presents a beautiful link between helically symmetric structure AE and toroidally symmetric structure GAM. Figure 7 shows that the frequency of AE decreases in time, caused by the density increase in accordance with Alfvén law f AE ∼ B t n −1/2 e . The circles indicate different branches of the Alfvén eigenmodes with different poloidal mode numbers m = 4 and m = 3, calculated from the data of a set of magnetic probes.
Symmetry 2021, 13, x FOR PEER REVIEW 12 of 49 Figure 8. Spectrogram of the coherency between the total HIBP current (plasma density) oscillations, measured in the plasma core at r/a = 0.16 and a MP signal [45]. Reproduced courtesy of IAEA. Copyright 2010 IAEA. Figure 9 shows an example of calculating the cross-phase spectrogram. We see that the magnitude of the cross-phase does not change, remaining equal to −π/2, although the mode frequency decreases twice with time, due to an increase in the density in this discharge. This constancy of the cross-phase values between the signals of the plasma density and the poloidal magnetic field indicates that the cross-phase is a stable individual characteristic for each Alfvén eigenmode. Figure 8. Spectrogram of the coherency between the total HIBP current (plasma density) oscillations, measured in the plasma core at r/a = 0.16 and a MP signal [45]. Reproduced courtesy of IAEA. Copyright 2010 IAEA. Figure 9 shows an example of calculating the cross-phase spectrogram. We see that the magnitude of the cross-phase does not change, remaining equal to −π/2, although the mode frequency decreases twice with time, due to an increase in the density in this discharge. This constancy of the cross-phase values between the signals of the plasma density and the poloidal magnetic field indicates that the cross-phase is a stable individual characteristic for each Alfvén eigenmode. Figure 9 shows an example of calculating the cross-phase spectrogram. We see that the magnitude of the cross-phase does not change, remaining equal to −π/2, although the mode frequency decreases twice with time, due to an increase in the density in this discharge. This constancy of the cross-phase values between the signals of the plasma density and the poloidal magnetic field indicates that the cross-phase is a stable individual characteristic for each Alfvén eigenmode. Figure 9. Spectrogram of the cross-phase between signals It (~ne) and poloidal magnetic field Bp (~ζd) measured in the same SV in discharge with AEs. Insets show the histograms of the cross-phase computed over the areas marked by rectangles; the insets show the histograms for the cross-phases [45]. Reproduced courtesy of IAEA. Copyright 2010 IAEA.

Bispectral Fourier Analysis of Plasma Oscillations
Bispectral analysis is widely used to study nonlinear three-wave interaction processes. The quadratic coefficient of cross-bicoherence for three independent quantities x(t), y(t) and z(t) is defined as follows: Figure 9. Spectrogram of the cross-phase between signals I t (~n e ) and poloidal magnetic field B p (~ζ d ) measured in the same SV in discharge with AEs. Insets show the histograms of the crossphase computed over the areas marked by rectangles; the insets show the histograms for the crossphases [45]. Reproduced courtesy of IAEA. Copyright 2010 IAEA.

Bispectral Fourier Analysis of Plasma Oscillations
Bispectral analysis is widely used to study nonlinear three-wave interaction processes. The quadratic coefficient of cross-bicoherence for three independent quantities x(t), y(t) and z(t) is defined as follows: where brackets denote the time averaging. Similar to the cross-phase, biphase is defined as an argument of the complex function: Auto-bicoherence for random variable x(t) is defined as: Statistically significant bicoherence b points to three-wave interaction between oscillations at frequencies f 1 , f 2 and f 3 = f 1 + f 2 . An example of calculating bicoherence and biphase is shown in Figure 10. The excess of the bicoherence coefficient over the background at the GAM frequency indicates a significant three-wave interaction between the GAM and the broadband turbulence of potential [46]. Note that, at the GAM frequencies the biphase has finite values, while outside them it has a stochastic character. Thus, analysis of HIBP data shows that GAM in T-10 linked with broadband turbulence [47,48], and is not excited by other mechanisms, such as, e.g., the effect of energetic particles.
lations at frequencies f1, f2 and f3 = f1 + f2. An example of calculating bicoherence and biphase is shown in Figure 10. The excess of the bicoherence coefficient over the background at the GAM frequency indicates a significant three-wave interaction between the GAM and the broadband turbulence of potential [46]. Note that, at the GAM frequencies the biphase has finite values, while outside them it has a stochastic character. Thus, analysis of HIBP data shows that GAM in T-10 linked with broadband turbulence [47,48], and is not excited by other mechanisms, such as, e.g., the effect of energetic particles.
(a) (b) Figure 10. Auto-bicoherency (a) and biphase (b) for the plasma potential in the ohmic discharge of the T-10 tokamak [46]. Reproduced courtesy of IAEA. Copyright 2017 IAEA.

Measurement of Turbulent Particle Flux on the TJ-II Stellarator
A fundamentally new element in the HIBP operation is the use of a multichannel energy analyzer, which allows one to measure both the potential and density in several sample volumes. Figure 11 shows a cross-sectional view of a multichannel analyzer installed in T-10. On the TJ-II stellarator, where two diagnostic systems HIBP1 and HIBP2 are currently operating (see Section 7), on the second HIBP2 system, the same analyzer is installed, somewhat reduced in size. Let us consider the operation of a multichannel analyzer using the simplest example of a two-channel analyzer with two slots and two detectors. Such an analyzer is installed on the HIBP1 complex of the TJ-II stellarator. The detector line for the TJ-II two-channel analyzer is shown in Figure 12.

Measurement of Turbulent Particle Flux on the TJ-II Stellarator
A fundamentally new element in the HIBP operation is the use of a multichannel energy analyzer, which allows one to measure both the potential and density in several sample volumes. Figure 11 shows a cross-sectional view of a multichannel analyzer installed in T-10. On the TJ-II stellarator, where two diagnostic systems HIBP1 and HIBP2 are currently operating (see Section 7), on the second HIBP2 system, the same analyzer is installed, somewhat reduced in size. Let us consider the operation of a multichannel analyzer using the simplest example of a two-channel analyzer with two slots and two detectors. Such an analyzer is installed on the HIBP1 complex of the TJ-II stellarator. The detector line for the TJ-II two-channel analyzer is shown in Figure 12.  The probing parameters are adjusted in such a way that both sample volumes lie on the same magnetic surface, and they are shifted poloidally by ∆x~1 cm, as shown in Figure 13. This method allows us to directly measure the poloidal component of the electric field E p by the potential difference, E p = (ϕ 1 − ϕ 2 )/∆x. The finite distance between the investigated volumes limits the poloidal wave vector of oscillations k p = k θ < 3 cm −1 . Direct measurement of E p allows us to define the radial velocity of E × B drift V r = E p /B t . Figure 11. Multi-slit energy analyzer for direct measurements of particle flux: 5-S-five entrance slits; D-detector; G-grid; GP-ground plate; HVP-high-voltage plate; W-adjustment window.

Figure 12.
Detector line for the two-slit energy analyzer [15]. The upper slit is marked in blue, the lower one-in red. Reproduced courtesy of IAEA. Copyright 2017 IAEA.
The probing parameters are adjusted in such a way that both sample volumes lie on the same magnetic surface, and they are shifted poloidally by ∆x ~1 cm, as shown in Figure  13. This method allows us to directly measure the poloidal component of the electric field Ep by the potential difference, Ep = (φ1 − φ2)/∆x. The finite distance between the investigated volumes limits the poloidal wave vector of oscillations kp = kθ < 3 cm −1 . Direct measurement of Ep allows us to define the radial velocity of E  B drift Vr = Ep/Bt.
As oscillations of the total particle current or plasma density are measured simultaneously with the potential oscillations, we may directly calculate the radial turbulent particle flux: Such use of HIBP on stellarators was first applied on the TJ-II. To find  As oscillations of the total particle current or plasma density are measured simultaneously with the potential oscillations, we may directly calculate the radial turbulent particle flux: Such use of HIBP on stellarators was first applied on the TJ-II. To find Γ E×B (t), the density fluctuations n e must be measured simultaneously with E p in the same spatial region. To analyze the frequency spectrum of the flux and its temporal dynamics, it is sufficient to measure the relative density fluctuations δn e (t) = I t (t)/I t . In the low-density case, to estimate the absolute values of the flux Γ E×B (t) we may use the measured value I t (t) as n e .
Symmetry 2021, 13, x FOR PEER REVIEW 15 of 49 Figure 13. Scheme of the poloidal electric field measurements in the plasma using the double-slit analyzer. The potential φ1 and φ2 are measured by the first and second slit correspondingly; ∆x is distance between the points of measurements; Ep is poloidal electric field; Bt is toroidal magnetic field, Vr is the velocity of radial drift in crossed Ep × Bt fields.
In the high-density case, the beam attenuation should be accounted.
, or the spectral function of flux is calculated as follows [49]: Here SnE is the Fourier cross-spectrum between density and poloidal electric field Ep oscillations. Figure 14 shows an example of the calculated turbulent particle flux in discharge with NBI heating. The flux, linked with broadband turbulence, demonstrates the intermittent nature. It consists of a stochastic sequence of controversially directed bursts, which are mostly directed outwards. The spectrograms show Alfvén eigenmodes in the form of a family of quasi-monochromatic oscillations with a frequency decreasing in time. We see Figure 13. Scheme of the poloidal electric field measurements in the plasma using the double-slit analyzer. The potential ϕ 1 and ϕ 2 are measured by the first and second slit correspondingly; ∆x is distance between the points of measurements; E p is poloidal electric field; B t is toroidal magnetic field, V r is the velocity of radial drift in crossed E p × B t fields.
In the high-density case, the beam attenuation should be accounted. The density should be normalized as n e = I tot /I tot · n e , where the oscillating density component I t /I t is measured with HIBP, and the normalized factor n e is measured by another diagnostic, e.g., by the interferometry or by Thomson scattering.
Fourier transform for the time-dependent flux Γ E×B (t), or the spectral function of flux is calculated as follows [49]: Here S nE is the Fourier cross-spectrum between density and poloidal electric field E p oscillations. Figure 14 shows an example of the calculated turbulent particle flux in discharge with NBI heating. The flux, linked with broadband turbulence, demonstrates the intermittent nature. It consists of a stochastic sequence of controversially directed bursts, which are mostly directed outwards. The spectrograms show Alfvén eigenmodes in the form of a family of quasi-monochromatic oscillations with a frequency decreasing in time. We see that AEs may contribute both to outward and inward flux, or also could produce no flux, depending on the phase relations between E p and the density oscillations [50]. Usually, AEs contribute to the turbulent particle flux Γ AE E×B , making up an appreciable fraction of the total flux Γ E×B . Figure 14 shows that Γ AE E×B can be much larger than the broadband turbulent flux Γ BB E×B in the same frequency range [51]. It was shown that the value of the turbulent flux is comparable to the value of the total particle flux; therefore, Γ E×B plays an essential role in the particle balance. It was shown that the transition to the improved confinement mode is characterized by suppression of broadband fluctuations of the plasma potential and density, as well as by a decrease Γ E×B both at the plasma edge and in the hot core [52].

Measurement of Turbulent Particle Flux in the T-10 Tokamak
In contrast to the TJ-II stellarator, the geometry of the T-10 machine gives more un certainties in determining the flux. The scheme of the experiment on measuring the flux with the 5-slit analyzer is shown in Figure 11. An example of calculating trajectories fo the standard discharge: Bt = 2.4 T, and energy Eb = 220 keV is shown in Figure 15. We can

Measurement of Turbulent Particle Flux in the T-10 Tokamak
In contrast to the TJ-II stellarator, the geometry of the T-10 machine gives more uncertainties in determining the flux. The scheme of the experiment on measuring the flux with the 5-slit analyzer is shown in Figure 11. An example of calculating trajectories for the standard discharge: B t = 2.4 T, and energy E b = 220 keV is shown in Figure 15. We can see that the potential difference ϕ k − ϕ j corresponds to the electric field E, directed along the segment k j . To measure the poloidal field, it is necessary to align the segment in the poloidal direction. This is achieved in the point of deepest penetration for the detector line of equal energy. [53]. Thus, we can measure the flux for each energy in one spatial point.   ϕ k is measured in the point S k by slit with number k = 1-5; E p is the poloidal component of the electric field required to calculate V r ; the radial E p × B t drift velocity in the crossed fields, and E r is its radial component. Figure 16 presents an example of the turbulent flux measurements in the standard T-10 discharge with ohmic and ECR heating. Only two frequency ranges with non-zero flux are seen in the spectral flux function: the quasi-coherent (QC) mode and the geodesic acoustic mode (GAM). In the remaining frequency ranges, the flux is stochastic. The QC mode dominates during the entire discharge, allowing particles to flow outward. The flux in the GAM frequency range Γ GAM E×B is by the factor of 5-6 less than the flux of the QC mode Γ QC E×B [54].
10 discharge with ohmic and ECR heating. Only two frequency ranges with non-zero flux are seen in the spectral flux function: the quasi-coherent (QC) mode and the geodesic acoustic mode (GAM). In the remaining frequency ranges, the flux is stochastic. The QC mode dominates during the entire discharge, allowing particles to flow outward. The flux in the GAM frequency range GAM E B  Γ is by the factor of 5-6 less than the flux of the QC mode QC E B  Γ [54].  Another example of the flux calculation for GAM is shown in Figure 17. Here we consider the ohmic T-10 discharge with slow density ramp-up. The density reaches the critical value at t = 600 ms, then the MHD instability develops, which leads to the disruption at t = 750 ms. Figure 17a shows that GAM dominates on the spectrogram of potential fluctuations, but its amplitude decreases with density increase. The spectrogram of fluctuations in the potential difference (Figure 17b) completely contrasts with Figure 17a, as in Figure 17b GAM does not dominate, but it is completely invisible due to the background noise. This observation agrees with the well-known fact that GAM on the potential has poloidal symmetry. Additionally, GAM was not observed on the spectrogram of flux [55].
An important conclusion follows from the above experiments: the radial turbulent particle flux at GAM frequencies is on the noise level or only slightly exceeds it due to the poloidal symmetry of potential oscillations at GAM frequencies.
in Figure 17b GAM does not dominate, but it is completely invisible due to the background noise. This observation agrees with the well-known fact that GAM on the potential has poloidal symmetry. Additionally, GAM was not observed on the spectrogram of flux [55].
An important conclusion follows from the above experiments: the radial turbulent particle flux at GAM frequencies is on the noise level or only slightly exceeds it due to the poloidal symmetry of potential oscillations at GAM frequencies.

Measurement of the GAM Poloidal Symmetry in the T-10 Tokamak
In the previous Subsection, we considered correlation measurements of potential, which are physically corrected as the potential measurements are completely local. In the low-density case, when the attenuation factors in the analysis of the variable component of the density may be neglected, the correlations of the density fluctuations may be also considered. The issue of the locality of density measurements is discussed in detail in the

Measurement of the GAM Poloidal Symmetry in the T-10 Tokamak
In the previous Subsection, we considered correlation measurements of potential, which are physically corrected as the potential measurements are completely local. In the low-density case, when the attenuation factors in the analysis of the variable component of the density may be neglected, the correlations of the density fluctuations may be also considered. The issue of the locality of density measurements is discussed in detail in the Reference [34]. Finally, as shown in Section 4, the density oscillations measurements are localized in the sample volume, so the density coherence is free from the non-local effects.
An example of a two-point measurement technique [56,57] is shown in Figure 18. The figure shows that the sample volumes are separated in the poloidal direction by 1.67-1.74 cm, but also have some radial shift of 0.81 cm.
The measurement results are shown in Figure 19. Time traces of the main discharge parameters are shown in Figure 19a. Figure 19b,c shows the low-frequency part of the power spectrograms of the potential oscillations measured using the center ϕ 3 (Figure 19b) and edge ϕ 5 (Figure 19c) slits. Figure 19d shows that the coherence coefficient between ϕ 3 and ϕ 5 in the GAM frequency domain reaches unity. Figure 19e shows that the phase shift between these poloidally shifted signals is zero, indicating that GAM is poloidal symmetric, i.e., m = 0.
Reference [34]. Finally, as shown in Section 4, the density oscillations measurements are localized in the sample volume, so the density coherence is free from the non-local effects.
An example of a two-point measurement technique [56,57] is shown in Figure 18. The figure shows that the sample volumes are separated in the poloidal direction by 1.671.74 cm, but also have some radial shift of 0.81 cm. Figure 18. Two-point technique of correlation measurements of the plasma electric potential using the 5-slit analyzer.
The measurement results are shown in Figure 19. Time traces of the main discharge parameters are shown in Figure 19a. Figure 19b,c shows the low-frequency part of the power spectrograms of the potential oscillations measured using the center φ3 ( Figure 19b) and edge φ5 (Figure 19c) slits. Figure 19d shows that the coherence coefficient between φ3 and φ5 in the GAM frequency domain reaches unity. Figure 19e shows that the phase shift between these poloidally shifted signals is zero, indicating that GAM is poloidal symmetric, i.e., m = 0.  Figure 20 shows the results of similar measurements carried out in the same mode on the radial interval of 18-26 cm. Figure shows that for closely located points 2 and 4, the cross-phase at the GAM frequency is close to zero over the entire measurement interval. The averaged error of cross-phase measurement, determined by the signal-to-noise ratio, is about 0.1 radians, increasing towards the edge to 0.2 radians. This gives the poloidal mode number m = 0 c with error ±1, and ±2 at the edge. Figure 21 shows the results of the   cross-phase at the GAM frequency is close to zero over the entire measurement interval. The averaged error of cross-phase measurement, determined by the signal-to-noise ratio, is about 0.1 radians, increasing towards the edge to 0.2 radians. This gives the poloidal mode number m = 0 c with error ±1, and ±2 at the edge. Figure 21 shows the results of the cross-phase measurements for the most distant points 1 and 5 in the same shot, but in a wider radial range. We see that the cross-phase between the signals from the 1st and 5th points at the GAM frequency is also close to zero over the entire measurement interval. The error of phase measurement is a constant of about 0.1 radians on average, and the error in determining the poloidal mode number m = 0 decreased due to an increase in the distance between points S 1 and S 5 on average to ±0.5, except for the edge, where it remained at a level of ±2 [58]. Thus, phase measurements of the potential confirmed the theoretically predicted poloidal symmetry of potential oscillations for GAM.

TM-4 HIBP
The first HIBP complex in Russia (USSR) was on the TM-4 tokamak, created early in the 1980s by Krupnik, Bondarenko and Nedzelskiy [59]. The scheme of the complex is shown in Figure 22.

TM-4 HIBP
The first HIBP complex in Russia (USSR) was on the TM-4 tokamak, created early in the 1980s by Krupnik, Bondarenko and Nedzelskiy [59]. The scheme of the complex is shown in Figure 22.

TM-4 HIBP
The first HIBP complex in Russia (USSR) was on the TM-4 tokamak, created early in the 1980s by Krupnik, Bondarenko and Nedzelskiy [59]. The scheme of the complex is shown in Figure 22.

TM-4 HIBP
The first HIBP complex in Russia (USSR) was on the TM-4 tokamak, created early in the 1980s by Krupnik, Bondarenko and Nedzelskiy [59]. The scheme of the complex is shown in Figure 22. The energy of secondary particles was determined by a multigrid ion analyzer. On TM-4, the first measurements of the potential profiles were carried out at a fixed position of the sample volume (SV), depending on the discharge parameters, such as the average chord density, plasma current and toroidal magnetic field. These measurements were carried out in a series of reproduced shots, in which the SV scans plasma from the center to the edge. It was shown that the plasma potential at the core has a negative sign and a scale about several hundred volts. For the first time, the link between the potential and the average plasma density was established: a higher density corresponds to a larger absolute value of negative potential. The measurements of the plasma density profile with a HIBP are described in Section 3.2. The TM-4 tokamak was terminated in 1982 due to the construction of the T-15 tokamak. Creation and scientific exploitation of TM-4 HIBP complex made the physical and technological basis for further development of the HIBP diagnostics towards the tool to study symmetric structures in toroidal plasmas.

T-10 HIBP
The HIBP complex on the T-10 tokamak was created in the mid-1980s by Krupnik, Bondarenko, Khrebtov and Nedzelskiy [18,60]. In contrast with the TM-4 complex, it used a Cs ion beam accelerated up to an energy of 160 keV, and a 'flat mirror' analyzer [61] with a higher energy resolution than the grid analyzer of TM-4. As the complex was implemented into the existing machine, there were no convenient ports to conduct the beam through the plasma. For injection we use a standard wide vertical port, and to pull out the secondary beam, we have to cut the vacuum chamber and weld an additional port to it [62]. Such a significant reconstruction of the machine became possible due to the successful HIBP operation on the TM-4 and the obtaining of important physical resultsfirst measurements of plasma potential profile. This important decision was adopted with support by Dnestrovskij, Razumova and Sokolov. The initial version of the HIBP complex is schematically shown in Figure 23.
secondary beam, we have to cut the vacuum chamber and weld an additional port to it [62]. Such a significant reconstruction of the machine became possible due to the successful HIBP operation on the TM-4 and the obtaining of important physical results-first measurements of plasma potential profile. This important decision was adopted with support by Dnestrovskij, Razumova and Sokolov. The initial version of the HIBP complex is schematically shown in Figure 23. Setting up the T-10 complex and conducting the ion beam into the analyzer encountered with an unexpected issue-a strong toroidal displacement of the beam, which hampered conducting the ion beam through the long narrow ports of the tokamak. To solve Setting up the T-10 complex and conducting the ion beam into the analyzer encountered with an unexpected issue-a strong toroidal displacement of the beam, which hampered conducting the ion beam through the long narrow ports of the tokamak. To solve it, for the first time in the HIBP operation, the toroidal correction of secondary particles was applied and a secondary ion beamline was installed outside the vacuum chamber. The passage of the beam into the analyzer through narrow diagnostic ports of the tokamak became possible thanks to the assistance of the T-10 leading experimentalist Dr. Bobrovskiy. At the initial stage of the HIBP operation, the measurements were carried out in a series of repeated tokamak shots, in which the SV position moved from the plasma edge towards the center. The limitations on the beam energy and ion mass allowed us to carry out measurements at a reduced magnetic field B t = 1.5 T in the range 0.6 < ρ SV < 1 [63][64][65].
Later on, the T-10 HIBP complex has permanently evolved in the following directions [66]:

1.
The Cs ions were replaced by heavier Tl ions. The energy range was consequently increased to 220, then to 280, 300, and finally to 330 kV. These changes allow to study discharges with microwave heating (ECRH) at 2.08 < B t < 2.5 T. In discharges with B t < 2.17 T, the radial region 0.2 < ρ SV < 1 was investigated; 2.
The ion beamlines were modernized in the following way: • Two serial pairs of steering plates were installed into the primary beamline. The first plates deflect the beam "counter-current" to the extreme toroidal position near the second plates. Thus, the beam approaches the second plates not axially, but maximally shifted. The second plates deflect the beam "co-current". As a result, a zigzag trajectory is formed, and it became possible to use the full toroidal width of the port, not half of it, as in the axial injection of the beam, and maximize the toroidal correction angle β. • The secondary beamline was significantly expanded. The toroidal correction plates were placed inside the vacuum chamber to increase the angle of toroidal correction of the secondary trajectories β 3 . Correcting plates were equipped with a continuous baking system that maintains their operating temperature at the level of 220-250 • C. This allows us to completely avoid the deposition of hydrocarbon films on the surface of the plates, and to exclude the appearance of high-voltage breakdowns between the plates or from the plate to the grounded T-10 wall. Immediately before the working shot of the tokamak, the baking was turned off and the baking circuit was disconnected in order to exclude the possibility of the interaction of the closed loop with the current and the confining magnetic fields of the machine. The electromagnetic forces resulting from this interaction could deform the plates, making it impossible for the secondary beam to pass into the analyzer. The described upgrades made it possible to expand the operational limits of measurements up to the operational limits of the tokamak 140 kA < I pl < 330 kA.

3.
The technique of scanning the entry angle of the particle beam into the plasma was implemented. Additionally, the vertical correction plates for secondary trajectories were installed to optimize the particle entry angle into the analyzer. It allowed us to obtain fragments of the radial profile per discharge; 4.
A control system for the beam was created, which allows us to select a complete set of control voltages in the primary and secondary ion beamlines in each subsequent tokamak shot based on the analysis of the beam position in the previous shot (injection and sweeping angles in primary and secondary beamlines α 1 , α 2 , α 3 , β 1 , β 2 , β 3 ). Its implementation radically simplified and accelerated the process of selecting control voltages, to significantly increase the accuracy of beam positioning, thereby providing systematic measurements both with a fixed SV position or with radial scan. Besides these two traditional HIBP operating modes, the developed system made it possible to operate in new non-standard modes. The most popular among them: • Periodic variation of the SV between two spatial positions ("colon"), • Periodic change of the SV between several positions ("multiple points"), • Alternation of the scan and at point measurements during one shot ("scan + point").

5.
A feedback for toroidal displacement was created, which adjust the toroidal correction voltage in the secondary beamline (toroidal angle β 3 ), depending on the measured toroidal displacement in the detector ζ d . This allows us to automate the selection of the toroidal correction voltage from shot to shot, as well as during one shot, especially in shots with changing plasma current (ramp-up, ramp-down) as well as the start of the plasma discharge with sharp increase of the plasma current. The feedback system for toroidal displacement is described in detail in the Reference [67]. 6.
The initial single-channel energy analyzer was first replaced by a two-channel and then a five-channel one. As a result, we can carry out correlation measurements over several spatial channels. In particular, in the region of maximal beam penetration into the plasma, it became possible to measure the poloidal potential correlations or E p , a turbulent particle flux Γ E×B , as well as the poloidal density correlations or rotation of turbulence. Closer to the plasma edge, it is possible to measure the velocity of radial correlations or the radial propagation of potential and density perturbations. 7.
An emitter-extractor unit has been developed, which makes it possible to significantly increase the primary beam current from 2-20 µA to 100-130 µA. Its use has expanded the allowable density limit for HIBP towards both high densities (plasma core) and also ultra-low densities at the plasma edge and in the scrape-off layer (SOL) [68].
The final version of the complex and its photo are shown in Figures 1 and 2. The T-10 tokamak was terminated in March 2018 due to the construction of the T-15MD tokamak. Creation and development of the T-10 HIBP complex provides comprehensive study of Geodesic Acoustic Mode (GAM)-the toroidally symmetric structures in tokamak plasmas, see Sections 5.1 and 6.3.

TJ-I HIBP
The HIBP complex for the TJ-I tokamak in the CIEMAT Science Center in Madrid, Spain, was created early in the 1990s by Krupnik, Bondarenko, Khrebtov and Nedzelskiy [69]. It operated with a Cs + ion beam, accelerated up to energy E b = 60 keV. The scheme of the complex is shown in Figure 24. Continuing the T-10 experiment, the TJ-I used the correcting toroidal plates for the secondary beam. As well as in T-10, the analyzer of the 'flat mirror' type was used on the TJ-I. The TJ-I tokamak had rectangular toroidal field coils. The TJ-I was the first, where the radial scanning method was used. The results of its application for E b = 20 keV are shown in Figure 25. We see that the plasma density has a bell-shaped profile with a pronounced maximum in the center, which shifts outward along the major radius during the discharge. The distribution of the plasma potential along the radius is inhomogeneous and also evolves with the density: the density decrease is accompanied by potential evolution towards the positive values.
The TJ-I tokamak was terminated in 1995 due to the construction of the TJ-II stellarator. Creation and scientific exploitation of TJ-I HIBP complex made the basis for the next-step TJ-II HIBP diagnostics, the most fruitful one in studies of the symmetric structures in toroidal plasmas.

TUMAN-3M HIBP
The HIBP complex on the small aspect ratio tokamak TUMAN-3M in the A.F. Ioffe Institute in St. Petersburg was created in the late 1990s by Krupnik, Bondarenko and Komarov [70]. It operates with a beam of K + ions accelerated up to the energy Eb = 60 keV. The scheme of the HIBP complex is shown in Figure 26. Figure shows that the detector line connects the center and edge of the plasma. It is suitable for determination of radial profiles over a total radial interval.

TUMAN-3M HIBP
The HIBP complex on the small aspect ratio tokamak TUMAN-3M in the A.F. Ioffe Institute in St. Petersburg was created in the late 1990s by Krupnik, Bondarenko and Komarov [70]. It operates with a beam of K + ions accelerated up to the energy Eb = 60 keV. The scheme of the HIBP complex is shown in Figure 26. Figure shows that the detector line connects the center and edge of the plasma. It is suitable for determination of radial profiles over a total radial interval.

TUMAN-3M HIBP
The HIBP complex on the small aspect ratio tokamak TUMAN-3M in the A.F. Ioffe Institute in St. Petersburg was created in the late 1990s by Krupnik, Bondarenko and Komarov [70]. It operates with a beam of K + ions accelerated up to the energy E b = 60 keV. The scheme of the HIBP complex is shown in Figure 26. Figure shows that the detector line connects the center and edge of the plasma. It is suitable for determination of radial profiles over a total radial interval. Figure 27 shows the equatorial projection of the calculated trajectories of the probing particles. We see that the plasma current shifts the trajectories, as well as that variation of the injection toroidal angle β effectively changes the displacement of the trajectories. Choosing β, we can place the detector line close to the meridian plane of the torus, z = 0. On the TUMAN-3M, it is possible to measure the evolution of the potential and plasma density at fixed spatial points.  Figure 27 shows the equatorial projection of the calculated trajectories of the probing particles. We see that the plasma current shifts the trajectories, as well as that variation of the injection toroidal angle β effectively changes the displacement of the trajectories. Choosing β, we can place the detector line close to the meridian plane of the torus, z = 0. On the TUMAN-3M, it is possible to measure the evolution of the potential and plasma density at fixed spatial points.

WEGA HIBP
Stellarator WEGA is a two-period torsatron with a low toroidal magnetic field and rotational transform ι/2π. The HIBP complex on the WEGA stellarator in the Max Planck Institute for Plasma Physics in Greifswald, Germany was created early in the 2000s by Krupnik,. It operates with a beam of Na + ions accelerated up to the energy E b = 50 keV. The scheme of the HIBP complex is shown in Figure 28, and the detector grid is shown in Figure 29.  Stellarator WEGA was terminated in 2006 due to the construction of the W7-X stellarator. The scientific exploitation of WEGA HIBP brings an additional experience of the HIBP operation in stellarators, so contributes to the next-step TJ-II HIBP diagnostics, the most fruitful one in studies of the symmetric structures in toroidal plasmas.

URAGAN-2M HIBP
The HIBP complex on URAGAN-2M stellarator at the Institute of Plasma Physics of the Kharkov Institute of Physics and Technology of the Academy of Sciences of Ukraine was created in the mid-2010th by Krupnik, Kozachek, Komarov, Zhezhera and Chmyga [74,75]. It operates with a Cs + ion beam accelerated to energies up to Eb = 70 keV. The scheme of the HIBP complex is shown in Figure 30a, and the equatorial projection of trajectories is shown in Figure 30b. URAGAN-2M was recommissioned in 2006 after a long  Stellarator WEGA was terminated in 2006 due to the construction of the W7-X stellarator. The scientific exploitation of WEGA HIBP brings an additional experience of the HIBP operation in stellarators, so contributes to the next-step TJ-II HIBP diagnostics, the most fruitful one in studies of the symmetric structures in toroidal plasmas.

URAGAN-2M HIBP
The HIBP complex on URAGAN-2M stellarator at the Institute of Plasma Physics of the Kharkov Institute of Physics and Technology of the Academy of Sciences of Ukraine was created in the mid-2010th by Krupnik, Kozachek, Komarov, Zhezhera and Chmyga [74,75]. It operates with a Cs + ion beam accelerated to energies up to Eb = 70 keV. The scheme of the HIBP complex is shown in Figure 30a, and the equatorial projection of tra- Stellarator WEGA was terminated in 2006 due to the construction of the W7-X stellarator. The scientific exploitation of WEGA HIBP brings an additional experience of the HIBP operation in stellarators, so contributes to the next-step TJ-II HIBP diagnostics, the most fruitful one in studies of the symmetric structures in toroidal plasmas.

URAGAN-2M HIBP
The HIBP complex on URAGAN-2M stellarator at the Institute of Plasma Physics of the Kharkov Institute of Physics and Technology of the Academy of Sciences of Ukraine was created in the mid-2010th by Krupnik, Kozachek, Komarov, Zhezhera and Chmyga [74,75]. It operates with a Cs + ion beam accelerated to energies up to E b = 70 keV. The scheme of the HIBP complex is shown in Figure 30a, and the equatorial projection of trajectories is shown in Figure 30b. URAGAN-2M was recommissioned in 2006 after a long pause caused by the reconstruction. It is currently undergoing the stage of adjustment and reaching the design parameters.

TJ-II HIBP-The Most Advanced Diagnostics to Study Plasma Symmetric Structures
The TJ-II stellarator-heliac has been successfully operating from the early 2000s to the present at the National Fusion Laboratory of the CIEMAT Science Center in Madrid, Spain. From the start of its operation till the start of the world largest stellarator W7-X in 2016, the TJ-II was the world's second largest stellarator after the LHD and the largest in Europe. Now it is the third stellarator after W7-X and LHD.
The initial version of the HIBP complex on the TJ-II was created in the early 2000s in collaboration with Krupnik, Nedzelskiy, Khrebtov and Bondarenko from the Kharkov Institute of Physics and Technology (Ukraine) and the National Research Center Kurchatov Institute (Russia) [76,77]. It operates with a Cs + ion beam accelerated up to the energy Eb = 150 keV.
Photos of the HIBP complex on the TJ-II stellarator are shown in Figure 31. Figure 32 presents the scheme of the HIBP complex, and the magnetic configuration with detector line, along which the measurements were performed, is shown in frame. Figure 30. Scheme of the HIBP on the URAGAN-2M stellarator. The optimal detector line is marked by blue asterisks, the green solid line is the trajectory; I is injector, D is detector, and P is plasma; (a) vertical plane, (b) equatorial plane.

TJ-II HIBP-The Most Advanced Diagnostics to Study Plasma Symmetric Structures
The TJ-II stellarator-heliac has been successfully operating from the early 2000s to the present at the National Fusion Laboratory of the CIEMAT Science Center in Madrid, Spain. From the start of its operation till the start of the world largest stellarator W7-X in 2016, the TJ-II was the world's second largest stellarator after the LHD and the largest in Europe. Now it is the third stellarator after W7-X and LHD.
The initial version of the HIBP complex on the TJ-II was created in the early 2000s in collaboration with Krupnik, Nedzelskiy, Khrebtov and Bondarenko from the Kharkov Institute of Physics and Technology (Ukraine) and the National Research Center Kurchatov Institute (Russia) [76,77]. It operates with a Cs + ion beam accelerated up to the energy E b = 150 keV.
Photos of the HIBP complex on the TJ-II stellarator are shown in Figure 31. Figure 32 presents the scheme of the HIBP complex, and the magnetic configuration with detector line, along which the measurements were performed, is shown in frame.
stitute of Physics and Technology (Ukraine) and the National Research Center Kurchatov Institute (Russia) [76,77]. It operates with a Cs + ion beam accelerated up to the energy Eb = 150 keV.
Photos of the HIBP complex on the TJ-II stellarator are shown in Figure 31. Figure 32 presents the scheme of the HIBP complex, and the magnetic configuration with detector line, along which the measurements were performed, is shown in frame.  The TJ-II stellarator turned out to be the unique machine, for which the HIBP complex was planned at the stage of machine design, rather than its operation, unlike other HIBPs. For this reason, HIBP design was optimized in advance aiming for the unique properties: the ability to obtain full profiles of plasma parameters by scanning voltage variation in one discharge.
An advanced HIBP in TJ-II have passed several steps to become a diagnostics capable to resolve the issue of the symmetry/asymmetry of the toroidal plasma parameters.
The first step was done in resolving the issue of LFS-HFS asymmetry. The detector line, started with radial interval −1 < ρ < 0 in the HFS and producing many important results in the characterization of the potential profile in biasing experiments [78,79], electron confinement study [80,81] and transition to the improved confinement regime [82,83], was then extended towards the LFS, 0 < ρ < 1, see the inset to Figure 32. This way the first experimental test for the potential profile symmetry was performed. The results are shown in Figure 33. This observation was the first experimental prove of the potential symmetry LFS-HFS in the toroidal devices. In this experiment, plasma density was varied by an order of magnitude, from 0.3 to 3 × 10 19 m −3 . It was shown that plasma potential evolves dramatically with the density raise, changing the shape and the sign [84]. In the whole observed domain of plasma parameters (density, temperature, collisionality, way of heating: ECRH, NBI) the evolving plasma potential is symmetric in terms of LFS-HFS. Later on, the symmetric structure of plasma potential profiles was proven for more extended domain of plasma parameters [85,86]. The TJ-II stellarator turned out to be the unique machine, for which the HIBP complex was planned at the stage of machine design, rather than its operation, unlike other HIBPs. For this reason, HIBP design was optimized in advance aiming for the unique properties: the ability to obtain full profiles of plasma parameters by scanning voltage variation in one discharge.
An advanced HIBP in TJ-II have passed several steps to become a diagnostics capable to resolve the issue of the symmetry/asymmetry of the toroidal plasma parameters.
The first step was done in resolving the issue of LFS-HFS asymmetry. The detector line, started with radial interval −1 < ρ < 0 in the HFS and producing many important results in the characterization of the potential profile in biasing experiments [78,79], electron confinement study [80,81] and transition to the improved confinement regime [82,83], was then extended towards the LFS, 0 < ρ < 1, see the inset to Figure 32. This way the first experimental test for the potential profile symmetry was performed. The results are shown in Figure 33. This observation was the first experimental prove of the potential symmetry LFS-HFS in the toroidal devices. In this experiment, plasma density was varied by an order of magnitude, from 0.3 to 3 × 10 19 m −3 . It was shown that plasma potential evolves dramatically with the density raise, changing the shape and the sign [84]. In the whole observed domain of plasma parameters (density, temperature, collisionality, way of heating: ECRH, NBI) the evolving plasma potential is symmetric in terms of LFS-HFS. Later on, the symmetric structure of plasma potential profiles was proven for more extended domain of plasma parameters [85,86]. Operating with the full detector line −1 < ρ < 1, on the TJ-II, the first measurements in stellarators were carried out, in which the dependence of the potential on the line-averaged density for a fixed position of the sample volume (SV) was studied, as well as the potential profiles and their dependence on heating methods and plasma parameters were measured [87,88]. At this stage there was made a substantial contribution to the study of helically symmetric structure-Alfven Eigenmodes. Figure 34 shows an example of the Helicity Induced Alfven Eigenmode (HAE) detected and identified by HIBP combined with the MHD-code AE3D [51]. The next principally important step in the widening of the diagnostic capabilities of HIBP was the measurement of the two-dimensional (2D) spatial distribution of the plasma potential and turbulence. This 2D map was obtained by the movement of the detector line across the plasma with the Eb variation [89,90].  Operating with the full detector line −1 < ρ < 1, on the TJ-II, the first measurements in stellarators were carried out, in which the dependence of the potential on the lineaveraged density for a fixed position of the sample volume (SV) was studied, as well as the potential profiles and their dependence on heating methods and plasma parameters were measured [87,88]. At this stage there was made a substantial contribution to the study of helically symmetric structure-Alfven Eigenmodes. Figure 34 shows an example of the Helicity Induced Alfven Eigenmode (HAE) detected and identified by HIBP combined with the MHD-code AE3D [51]. The next principally important step in the widening of the diagnostic capabilities of HIBP was the measurement of the two-dimensional (2D) spatial distribution of the plasma potential and turbulence. This 2D map was obtained by the movement of the detector line across the plasma with the E b variation [89,90]. Operating with the full detector line −1 < ρ < 1, on the TJ-II, the first measurements in stellarators were carried out, in which the dependence of the potential on the line-averaged density for a fixed position of the sample volume (SV) was studied, as well as the potential profiles and their dependence on heating methods and plasma parameters were measured [87,88]. At this stage there was made a substantial contribution to the study of helically symmetric structure-Alfven Eigenmodes. Figure 34 shows an example of the Helicity Induced Alfven Eigenmode (HAE) detected and identified by HIBP combined with the MHD-code AE3D [51]. The next principally important step in the widening of the diagnostic capabilities of HIBP was the measurement of the two-dimensional (2D) spatial distribution of the plasma potential and turbulence. This 2D map was obtained by the movement of the detector line across the plasma with the Eb variation [89,90].  To construct a 2D map, the systematic measurements were performed in a series of the reproducible shots. The results are presented in Figure 35. This 2D map of plasma potential is the first experimental prove of the potential poloidal symmetry in the toroidal devices. Figure 35a shows that the contours of the plasma potential (equipotentials) are well in line with the magnetic flux surfaces [91]. In contrast, the 2D map of the plasma potential and density turbulence is not poloidally symmetric, it has higher level at the HFS, as presented in Figure 35b. To construct a 2D map, the systematic measurements were performed in a series of the reproducible shots. The results are presented in Figure 35. This 2D map of plasma potential is the first experimental prove of the potential poloidal symmetry in the toroidal devices. Figure 35a shows that the contours of the plasma potential (equipotentials) are well in line with the magnetic flux surfaces [91]. In contrast, the 2D map of the plasma potential and density turbulence is not poloidally symmetric, it has higher level at the HFS, as presented in Figure 35b. It is worthwhile to highlight that these observations are contrasting to the gyrokinetic expectations for toroidal devices, always predicting the contrary: higher level of turbulence at the LFS [92]. This experimental test should help the theorists to look deeper at the theoretical basis for the modelling.
The final step of the paramount importance is a duplication of the HIBP diagnostics in one machine and creation of the Dual HIBP [93,94]. This system is aimed for the study of the long-range toroidal correlations in plasma parameters, and, more specifically, zonal flow [95,96]. It consists of the former complex, now called HIBP-1 and the new HIBP-2, separated toroidally in a quarter of the torus as presented in Figure 36. The photo of the Dual HIBP is shown in Figure 37. It is worthwhile to highlight that these observations are contrasting to the gyrokinetic expectations for toroidal devices, always predicting the contrary: higher level of turbulence at the LFS [92]. This experimental test should help the theorists to look deeper at the theoretical basis for the modelling.
The final step of the paramount importance is a duplication of the HIBP diagnostics in one machine and creation of the Dual HIBP [93,94]. This system is aimed for the study of the long-range toroidal correlations in plasma parameters, and, more specifically, zonal flow [95,96]. It consists of the former complex, now called HIBP-1 and the new HIBP-2, separated toroidally in a quarter of the torus as presented in Figure 36. The photo of the Dual HIBP is shown in Figure 37.  Long-range toroidal correlations were found with the Dual HIBP in low-frequency plasma potential oscillations. Their radial distribution is presented in Figure 38. Figure 38 shows that there is a statistically significant coherence between plasma potential oscillations, measured by HIBP1 and HIBP2 [97,98]. The only low-frequency (<30 kHz) oscillations are coherent and having a zero cross-phase. Therefore, the observed phenomenon presents toroidally symmetric (n = 0) low-frequency potential perturbations, which is, by definition, the zonal flow [99].  Long-range toroidal correlations were found with the Dual HIBP in low-frequency plasma potential oscillations. Their radial distribution is presented in Figure 38. Figure 38 shows that there is a statistically significant coherence between plasma potential oscillations, measured by HIBP1 and HIBP2 [97,98]. The only low-frequency (<30 kHz) oscillations are coherent and having a zero cross-phase. Therefore, the observed phenomenon presents toroidally symmetric (n = 0) low-frequency potential perturbations, which is, by definition, the zonal flow [99]. Long-range toroidal correlations were found with the Dual HIBP in low-frequency plasma potential oscillations. Their radial distribution is presented in Figure 38. Figure 38 shows that there is a statistically significant coherence between plasma potential oscillations, measured by HIBP1 and HIBP2 [97,98]. The only low-frequency (<30 kHz) oscillations are coherent and having a zero cross-phase. Therefore, the observed phenomenon presents toroidally symmetric (n = 0) low-frequency potential perturbations, which is, by definition, the zonal flow [99]. In addition to study zonal flows [100], the double HIBP was also focused on the longrange correlations driven by helically symmetric structures-Alfvén eigenmodes [101][102][103][104][105] and induced by pellet injection [106,107].
The development and scientific exploitation of TJ-II dual HIBP have shown an importance of the symmetry and symmetry braking studies at the flux surface in toroidal devices. The results obtained erects the interest to widen such studies in other tokamaks and stellarators. The next section describes the HIBP projects under consideration for various devices in operation or in construction.

Future Prospects to Study of Symmetric Structures in Toroidal Plasmas-Conceptual Design of the HIBP Diagnostics for Various Toroidal Devices
The technical demands to the future HIBP projects are basically formulated in Sections 2-4. An experience of the operation and physical results obtained during the HIBP operation, discussed in Sections 5-7, form the new and more challenging level for the new projects. It consists of the following items: a) Simultaneous measurements of the all three signals for potential, density, magnetic oscillations for comprehensive analysis of the plasma phenomena, including turbulence; b) Multichannel measurements for correlation studies, including plasma turbulence poloidal rotation and radial propagation, and turbulent particle flux; c) Maximal extension of the detector grid for 2D mapping and study of the poloidal symmetry/symmetry braking; d) Creation of the dual HIBP system for study of toroidal/helical symmetry.
The development and scientific exploitation of TJ-II dual HIBP have shown an importance of the symmetry and symmetry braking studies at the flux surface in toroidal devices. The results obtained erects the interest to widen such studies in other tokamaks and stellarators. The next section describes the HIBP projects under consideration for various devices in operation or in construction.

Future Prospects to Study of Symmetric Structures in Toroidal Plasmas-Conceptual Design of the HIBP Diagnostics for Various Toroidal Devices
The technical demands to the future HIBP projects are basically formulated in Sections 2-4. An experience of the operation and physical results obtained during the HIBP operation, discussed in Sections 5-7, form the new and more challenging level for the new projects. It consists of the following items: (a) Simultaneous measurements of the all three signals for potential, density, magnetic oscillations for comprehensive analysis of the plasma phenomena, including turbulence; (b) Multichannel measurements for correlation studies, including plasma turbulence poloidal rotation and radial propagation, and turbulent particle flux; (c) Maximal extension of the detector grid for 2D mapping and study of the poloidal symmetry/symmetry braking; (d) Creation of the dual HIBP system for study of toroidal/helical symmetry.
The projects briefly presented below are oriented to address the critical issues of HIBP compatibility to the machines, and HIBP capability to reply to the new challenges (a)-(d).

HIBP Design for the TCABR Tokamak
The conceptual design of the HIBP complex for TCABR tokamak at the Plasma Physics Laboratory of the University of São Paulo in Brazil was developed in the mid-2000s by Kuznetsov and Krupnik on the initiative of Severo and Nasimento [108]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated up to the energy E b = 105 keV. The TCABR tokamak has rectangular toroidal field coils. The scheme of the HIBP complex is shown in Figure 39a, and the detector grid is shown in Figure 39b. At the TCABR tokamak, the polarization of the plasma edge (biasing) and heating by Alfvén waves [109] are studied. For both experimental programs, 2D data on the plasma potential and its fluctuations are required.

HIBP Design for the TCABR Tokamak
The conceptual design of the HIBP complex for TCABR tokamak at the Plasma Physics Laboratory of the University of São Paulo in Brazil was developed in the mid-2000s by Kuznetsov and Krupnik on the initiative of Severo and Nasimento [108]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated up to the energy Eb = 105 keV. The TCABR tokamak has rectangular toroidal field coils. The scheme of the HIBP complex is shown in Figure 39a, and the detector grid is shown in Figure 39b. At the TCABR tokamak, the polarization of the plasma edge (biasing) and heating by Alfvén waves [109] are studied. For both experimental programs, 2D data on the plasma potential and its fluctuations are required.

HIBP Design for the Globus-M2 Tokamak
The conceptual design of the HIBP complex for the Globus-M2 spherical tokamak with an extremely small aspect ratio (R/a = 0.36 m/0.24 m = 1.5), operating in the A.F. Ioffe Physical Institute, St. Petersburg, was based at the initial design studies for the earlier version Globus-M [70,110]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated up to Eb ≤ 45 keV. Figure 40 shows the scheme of the HIBP complex. Figure 40a shows both the vertical and the equatorial projections of the calculated trajectories of the probing particles. Figure 40b shows the detector grid aiming for the study of the 2D structures at the outer half of plasma radius in the LFS.

HIBP Design for the Globus-M2 Tokamak
The conceptual design of the HIBP complex for the Globus-M2 spherical tokamak with an extremely small aspect ratio (R/a = 0.36 m/0.24 m = 1.5), operating in the A.F. Ioffe Physical Institute, St. Petersburg, was based at the initial design studies for the earlier version Globus-M [70,110]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated up to E b ≤ 45 keV. Figure 40 shows the scheme of the HIBP complex. Figure 40a shows both the vertical and the equatorial projections of the calculated trajectories of the probing particles. Figure 40b shows the detector grid aiming for the study of the 2D structures at the outer half of plasma radius in the LFS. At the Globus-M2 tokamak, studies of plasma confinement under conditions of an extremely small aspect ratio and plasma heating with NBI and lower hybrid waves are conducted with focus to Alfvén eigenmodes. Most of experimental programs require data on the plasma potential and turbulence.

HIBP Design for the COMPASS Tokamak
The small aspect ratio COMPASS tokamak, operating at the Institute of Plasma Physics of the Czech Academy of Sciences in Prague, is used to study the L-H transition in discharges with ohmic and NBI plasma heating. The studies are focused to Alfvén eigenmodes and GAMs [111,112]. These experimental programs require data on the plasma potential and its fluctuations.
Conceptual design of the HIBP for the COMPASS tokamak (R = 0.56 m, a = 0.23 m, Ip < 360 kA, Bt ≤ 2.1, PNBI (at 40 keV) = 2 × 0.3 MW) with a, was created on the initiative of Prof. J. Stökel. COMPASS has a single-null divertor plasma configuration with elongation k = 1.8, triangularity ∆ = 0.2, and horizontal plasma size of 16 cm. The detector line, calculated in the classical probing scheme for and rectangular toroidal field coils, does not pass through the plasma center. Therefore, we consider another possible combination of probing ports, with the injection of primary particles through the lower inclined port and registration of secondary particles through the upper inclined port. The calculation and design results are shown in Figure 41. We see that the detector line connects the center and the plasma edge. Thus, this probing scheme is preferable. However, for realization of this scheme, the direction of the magnetic field Bt should be reversed. At the Globus-M2 tokamak, studies of plasma confinement under conditions of an extremely small aspect ratio and plasma heating with NBI and lower hybrid waves are conducted with focus to Alfvén eigenmodes. Most of experimental programs require data on the plasma potential and turbulence.

HIBP Design for the COMPASS Tokamak
The small aspect ratio COMPASS tokamak, operating at the Institute of Plasma Physics of the Czech Academy of Sciences in Prague, is used to study the L-H transition in discharges with ohmic and NBI plasma heating. The studies are focused to Alfvén eigenmodes and GAMs [111,112]. These experimental programs require data on the plasma potential and its fluctuations.
Conceptual design of the HIBP for the COMPASS tokamak (R = 0.56 m, a = 0.23 m, I p < 360 kA, B t ≤ 2.1, P NBI (at 40 keV) = 2 × 0.3 MW) with a, was created on the initiative of Prof. J. Stökel. COMPASS has a single-null divertor plasma configuration with elongation k = 1.8, triangularity ∆ = 0.2, and horizontal plasma size of 16 cm. The detector line, calculated in the classical probing scheme for and rectangular toroidal field coils, does not pass through the plasma center. Therefore, we consider another possible combination of probing ports, with the injection of primary particles through the lower inclined port and registration of secondary particles through the upper inclined port. The calculation and design results are shown in Figure 41. We see that the detector line connects the center and the plasma edge. Thus, this probing scheme is preferable. However, for realization of this scheme, the direction of the magnetic field B t should be reversed. Symmetry 2021, 13, x FOR PEER REVIEW 37 of 49 Figure 41. Design of the heavy ion beam probing scheme for the COMPASS tokamak. Figure 42 shows a detector grid corresponding to the scheme in Figure 41.

HIBP Design for the TCV Tokamak
The TCV (Tokamak á Configuracion Variable) tokamak at the Plasma Physics Laboratory of the Federal Polytechnic School in Lausanne, Switzerland, is used to study dependencies of plasma confinement depending on the shape of the magnetic configuration, the L-H transition under conditions of ohmic and NBI plasma heating. Particular attention is paid to the study of Alfvén eigenmodes and GAM. These experimental programs require data on the plasma potential and its fluctuations.
The conceptual design of the HIBP complex for the TCV was developed on the initiative of Tonetti and Weisen [113]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated to Eb = 105 keV. The TCV has rectangular toroidal field coils. The

HIBP Design for the TCV Tokamak
The TCV (Tokamak á Configuracion Variable) tokamak at the Plasma Physics Laboratory of the Federal Polytechnic School in Lausanne, Switzerland, is used to study dependencies of plasma confinement depending on the shape of the magnetic configuration, the L-H transition under conditions of ohmic and NBI plasma heating. Particular attention is paid to the study of Alfvén eigenmodes and GAM. These experimental programs require data on the plasma potential and its fluctuations.
The conceptual design of the HIBP complex for the TCV was developed on the initiative of Tonetti and Weisen [113]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated to Eb = 105 keV. The TCV has rectangular toroidal field coils. The

HIBP Design for the TCV Tokamak
The TCV (Tokamak á Configuracion Variable) tokamak at the Plasma Physics Laboratory of the Federal Polytechnic School in Lausanne, Switzerland, is used to study dependencies of plasma confinement depending on the shape of the magnetic configuration, the L-H transition under conditions of ohmic and NBI plasma heating. Particular attention is paid to the study of Alfvén eigenmodes and GAM. These experimental programs require data on the plasma potential and its fluctuations.
The conceptual design of the HIBP complex for the TCV was developed on the initiative of Tonetti and Weisen [113]. It is assumed that the diagnostics will operate with a Tl ion beam accelerated to E b = 105 keV. The TCV has rectangular toroidal field coils. The scheme of the HIBP complex is shown in Figure 43a, and the detector grid is shown in Figure 43b. scheme of the HIBP complex is shown in Figure 43a, and the detector grid is shown in Figure 43b.

HIBP Design for the MAST Tokamak
The MAST (Mega-Ampere Spherical Torus) tokamak at the Culham Science Center in Culham, UK, is used to study the plasma confinement in the configuration of a spherical tokamak, L-H transition under conditions of ohmic and NBI plasma heating, and divertor physics. Particular attention is paid to the study of the Alfvén eigenmodes and Edge Localized Modes (ELMs). These experimental programs require data on the plasma potential and its fluctuations.
The conceptual design of the HIBP complex for MAST spherical tokamak was developed in on the initiative of Winsor and Sharapov [114]. It is assumed that the diagnostics will work with the Tl ion beam accelerated to Eb = 105 keV. The MAST tokamak has parameters R = 0.85 m, a = 0.65 m, Bt = 0.6 T, Ipl = 1 MA, an extremely small aspect ratio and rectangular toroidal field coils. The vertical section of the tokamak vacuum chamber is also rectangular. The scheme of HIBP is shown in Figure 44a [115], and one of the probing options and the detector grid is shown in Figure 44b. Note that several trajectories (blue) have intersections that complicate interpretation of HIBP measurements, while other trajectories (red) have not intersections. The similar situation was in COMPASS (Figure 41).

HIBP Design for the MAST Tokamak
The MAST (Mega-Ampere Spherical Torus) tokamak at the Culham Science Center in Culham, UK, is used to study the plasma confinement in the configuration of a spherical tokamak, L-H transition under conditions of ohmic and NBI plasma heating, and divertor physics. Particular attention is paid to the study of the Alfvén eigenmodes and Edge Localized Modes (ELMs). These experimental programs require data on the plasma potential and its fluctuations.
The conceptual design of the HIBP complex for MAST spherical tokamak was developed in on the initiative of Winsor and Sharapov [114]. It is assumed that the diagnostics will work with the Tl ion beam accelerated to E b = 105 keV. The MAST tokamak has parameters R = 0.85 m, a = 0.65 m, B t = 0.6 T, I pl = 1 MA, an extremely small aspect ratio and rectangular toroidal field coils. The vertical section of the tokamak vacuum chamber is also rectangular. The scheme of HIBP is shown in Figure 44a [115], and one of the probing options and the detector grid is shown in Figure 44b

HIBP design for the T-15MD tokamak
The T-15 MD tokamak in the National Research Center Kurchatov Institute, Moscow, Russia, is commissioned in 2021. T-15 MD intended to study the plasma confinement in conditions of a strong magnetic field (Bt  2 T) combined with a small aspect ratio (R/a = 1.48 m/0.67 m). It is planned to use ohmic, ECRH, ICRH and NBI plasma heating and current drive (CD). Particular attention will be paid to the study of the role of the electric field in L-H transitions, Alfvén eigenmodes, and GAMs. These experimental programs require 2D data on the plasma potential, density, magnetic field and its fluctuations.
Currently, conceptual design of the HIBP complex for T-15 MD is developed as a part of the diagnostic complex of the machine [116]. The HIBP will use the Tl ion beam with Eb = 400 keV. T-15 MD magnetic system was taken into account in the detailed trajectory calculations [117]. The scheme of the HIBP complex is shown in Figure 45a [118], and the detector grid is shown in Figure 45b.

HIBP Design for the T-15MD Tokamak
The T-15 MD tokamak in the National Research Center Kurchatov Institute, Moscow, Russia, is commissioned in 2021. T-15 MD intended to study the plasma confinement in conditions of a strong magnetic field (B t ≤ 2 T) combined with a small aspect ratio (R/a = 1.48 m/0.67 m). It is planned to use ohmic, ECRH, ICRH and NBI plasma heating and current drive (CD). Particular attention will be paid to the study of the role of the electric field in L-H transitions, Alfvén eigenmodes, and GAMs. These experimental programs require 2D data on the plasma potential, density, magnetic field and its fluctuations.
Currently, conceptual design of the HIBP complex for T-15 MD is developed as a part of the diagnostic complex of the machine [116]. The HIBP will use the Tl ion beam with E b = 400 keV. T-15 MD magnetic system was taken into account in the detailed trajectory calculations [117]. The scheme of the HIBP complex is shown in Figure 45a [118], and the detector grid is shown in Figure 45b.
In the future, it is planned to upgrade diagnostics and create a double HIBP, consisting of two identical complexes, shifted along the torus by about 90 • . Thus, the dual HIBP on the T-15 MD tokamak will be basically similar to the dual HIBP on the TJ-II stellarator. The main task of the dual HIBP on T-15 MD is to study zonal flows and long-range correlations of potential and density at the frequencies of GAMs [119,120] and Alfvén eigenmodes [121] including the effects of the ECRH and CD [122]. On top of that, the special attention will be dedicated to the ICRH and CD [123,124]. In the future, it is planned to upgrade diagnostics and create a double HIBP, consisting of two identical complexes, shifted along the torus by about 90. Thus, the dual HIBP on the T-15 MD tokamak will be basically similar to the dual HIBP on the TJ-II stellarator. The main task of the dual HIBP on T-15 MD is to study zonal flows and long-range correlations of potential and density at the frequencies of GAMs [119,120] and Alfvén eigenmodes [121] including the effects of the ECRH and CD [122]. On top of that, the special attention will be dedicated to the ICRH and CD [123,124].

HIBP Design for the W7-X Stellarator
The conceptual design of the HIBP complex for the world's largest stellarator W7-X was developed in the early 2000s. This installation has the following design parameters: Bt = 3 T, R = 5.5 m, a = 0.53 m. W7-X was commissioned at the end of 2015 in the Max Planck Institute for Plasma Physics in Greifswald, Germany. W7-X intended to study the plasma confinement in the optimized 3D magnetic configuration. It uses microwave and NBI for plasma heating. Particular attention is focused on study the role of the electric field at L-H transitions and Alfvén eigenmodes. These experimental programs require data on the plasma potential and its fluctuations. A preliminary calculation of trajectories for the proposed W7-X plasma probing scheme is shown in Figure 46 [125,126].

HIBP Design for the W7-X Stellarator
The conceptual design of the HIBP complex for the world's largest stellarator W7-X was developed in the early 2000s. This installation has the following design parameters: B t = 3 T, R = 5.5 m, a = 0.53 m. W7-X was commissioned at the end of 2015 in the Max Planck Institute for Plasma Physics in Greifswald, Germany. W7-X intended to study the plasma confinement in the optimized 3D magnetic configuration. It uses microwave and NBI for plasma heating. Particular attention is focused on study the role of the electric field at L-H transitions and Alfvén eigenmodes. These experimental programs require data on the plasma potential and its fluctuations. A preliminary calculation of trajectories for the proposed W7-X plasma probing scheme is shown in Figure 46 [125,126]. (2)-fan of secondary trajectories directed to the analyzer through the port N11; 3-total fan of secondary trajectories.

HIBP Design for the International Experimental Tokamak Reactor ITER
The ITER, which is currently being assembled in the city of Cadarache in the south of France, aims to achieve the fusion power amplification factor Q = 5. The ITER diagnostics are subordinated to this engineering task. They are primarily focused to achieving the design parameters of the machine and provide the safety of its operation. So, HIBP is not

HIBP Design for the International Experimental Tokamak Reactor ITER
The ITER, which is currently being assembled in the city of Cadarache in the south of France, aims to achieve the fusion power amplification factor Q = 5. The ITER diagnostics are subordinated to this engineering task. They are primarily focused to achieving the design parameters of the machine and provide the safety of its operation. So, HIBP is not included in the list of top-priority diagnostics for ITER, but it may be used in studying the physics of the L-H transition. Therefore, the HIBP allows one to study the periphery of the ITER plasma. Estimations have shown that HIBP will operate with a Tl ion beam accelerated to the energy of E b = 3.5 MeV [127]. This is a standard energy range for largescale machines in operation [19]. ITER has ovoid-shaped toroidal field coils. The scheme of the HIBP complex is shown in Figure 47, and the detector grid is shown in Figure 48a. In the standard ITER mode, the nominal plasma current will lead to a noticeable toroidal shift of the trajectories of the probing particles, shown in Figure 48b. We see that with the initial toroidal shift z I = 0 at the injection point I and the value of the toroidal angle β = 7 • , the toroidal coordinate of the detection point D, z D = 10 cm, and easily fits into the dimensions of the horizontal port. Because in ITER the discharge will transit to the H-mode that is always accompanied by change in the electric field [52,128,129], the measurements of the electric potential and turbulence will be an important issue for the physics of ITER. Recently, it was first shown by HIBP in the T-10 tokamak that in the discharges with low collisionality-which mimic the ITER plasma condition-the core plasma potential has positive values and an enhanced level of turbulence, so the prediction for the ITER plasma was made [130]. Verification of this prediction could be a very important and challenging task for physics studies in ITER.

Summary
Over more than 30 years of development, diagnostics of fusion plasma using a beam of heavy ions has become a powerful and multifunctional tool for physical research at the magnetic confinement devices. Starting with a time-averaged potential signal measured in the single spatial point, it becomes the diagnostic for studying 2D spatial distributions of plasma potential and turbulence, and even 3D characteristics of turbulence. On modern operating medium-scale machines, this diagnostic makes a significant contribution to reveal the role of electric fields and turbulence in confinement and the transport of plasma energy and particles, physics of fast particles and Alfvén eigenmodes. These advances have allowed HIBP to study the spatial and frequency characteristics of toroidally symmetric structures: geodesic acoustic modes and zonal flows, and helically symmetric

Summary
Over more than 30 years of development, diagnostics of fusion plasma using a beam of heavy ions has become a powerful and multifunctional tool for physical research at the magnetic confinement devices. Starting with a time-averaged potential signal measured in the single spatial point, it becomes the diagnostic for studying 2D spatial distributions of plasma potential and turbulence, and even 3D characteristics of turbulence. On modern operating medium-scale machines, this diagnostic makes a significant contribution to reveal the role of electric fields and turbulence in confinement and the transport of plasma energy and particles, physics of fast particles and Alfvén eigenmodes. These advances have allowed HIBP to study the spatial and frequency characteristics of toroidally symmetric structures: geodesic acoustic modes and zonal flows, and helically symmetric structures: Alfvén Eigenmodes. The future studies will be focused on the further development of operating HIBP capabilities and the installations to the other considered toroidal devices.