A model of non-Maxwellian electron distribution function for the analysis of ECE data in JET discharges

. Recent experiments performed in JET at high level of plasma heating, in preparation of, and during the DT campaign have shown significant discrepancies between electron temperature measurements by Thomson Scattering (TS) and Electron Cyclotron Emission (ECE). In order to perform a systematic analysis of this phenomenon, a simple model of bipolar distortion of the electron distribution function has been developed, allowing analytic calculation of the EC emission and absorption coefficients. Extensive comparisons of the modelled ECE spectra (at both the 2 nd and the 3 rd harmonic extraordinary mode) with experimental measurements display good agreement when bulk electron distribution distortions around 1-2 times the electron thermal velocity are used and prove useful for a first level of analysis of this effect.


Introduction
Discrepancies between electron temperature measurements by Thomson Scattering (TS) and Electron Cyclotron Emission (ECE) have been often observed in high-temperature tokamak plasmas, in particular on TFTR [1], JET [2] and FTU [3]. Such observations, made on different machines, by different types of instruments, using different calibration methods, are too ubiquitous to be ascribed to instrumental effects; they rather call for explanations based on physics phenomena. The hypothesis that the discrepancy could be associated to non-Maxwellian bulk electron distributions has been put forward in the past [4,5] and appears as a plausible explanation in the case of a plasma heated by EC waves, as FTU [3]. For TFTR and JET, electron heating rather takes place because of the interaction of the electron distribution either with a fast ion tail driven by Neutral Beam Injection (NBI) and/or Ion Cyclotron Resonance Heating (ICRH), or by energetic alpha particles produced by fusion reactions in DT (Deuterium-Tritium) plasmas. Two mechanisms are known to produce small bipolar distortions of the electron distribution in the presence of energetic ions: collisional * Corresponding author: gerardo.giruzzi@cea.fr relaxation [6] or Landau damping of kinetic Alfvén waves, as observed in the magnetosheath [7,8].
Recent experiments performed at JET at high level of plasma heating, in preparation of, and during the DT campaign have shown again TS-ECE discrepancies on an extensive database [9]. ECE is observed to be higher or lower than TS, depending on the plasma scenario. Moreover, ECE measured by a Martin-Puplett interferometer on a broad frequency range displays differences between 2 nd and 3 rd harmonics extraordinary (X) mode, which, at high temperatures (> 4 keV) and high densities are expected to yield the same radiation temperature. In order to perform a systematic analysis of this effect, a simple model of bipolar distortion of the electron distribution function has been developed, allowing analytic calculation of the EC emission and absorption coefficients. Bulk electron distribution distortions around 1-2 times the electron thermal velocity are considered for a first level of analysis of this effect. In this paper, the model is described in Sec. 2. Comparisons of the modelled ECE spectra (at both the 2 nd and the 3 rd harmonic) with experimental measurements are presented in Sec. 3. Conclusions and possible impact on ITER are presented in Sec. 4.

Electron distribution function model
A toy model of isotropic perturbation f1(p) of the electron distribution function f(p), where p is the modulus of the electron momentum, is developed as follows. We take the relativistic Maxwellian as the unperturbed distribution: m is the electron rest mass, c the speed of light, Te the electron temperature and K2 the modified Bessel function of the second kind. The perturbed distribution function is defined as: and a suitable bipolar isotropic form of the perturbation f1 is given as a function of three parameters f0, p0 and  by: Various types of anisotropic forms can also be defined by multiplying f1u by functions of the pitchangle , e.g., Two examples of isotropic and anisotropic perturbed distribution functions of this kind are shown in Fig. 1. All these functions allow analytical calculations of the electron cyclotron absorption () and emission () coefficients for perpendicular propagation. The general expressions of the absorption and emission coefficients for an arbitrary electron distribution function, as momentum-space integrals, are well known and available in a number of papers. Starting, for instance, from Eqs. 10-13 of Ref. [10], for the extraordinary wave and n > 1, these expressions are given by: where n is the harmonic number,  the wave frequency (times 2), p and c the plasma and the electron cyclotron frequencies, Nx the cold refractive index and  the cold plasma dielectric tensor; (x,y) is the Beta function and H the Heaviside function.
The function Qn depends on the choice of the perturbed distribution function as follows: In the following, only the isotropic perturbation f1u is considered.
Using these expressions of the emission and absorption coefficients, the radiation temperature, measured along a line of sight in the equatorial plane (as it is approximately the case for JET), is given by: where R is the major radius coordinate and a, R0 are the minor and major radii. An example of the absorption and emission coefficients, as well as of the radiation function = (− ∫ ) , i.e., the integrand appearing in the radiation temperature expression, at various wave frequencies and for typical JET parameters is presented in Fig. 2. This figure shows that  and  are broad functions of R, nevertheless, emission is well localised in space owing to the exponential reabsorption term that multiplies  in the radiation function expression.
Owing to the resonance condition  = nc(R)/, for a given wave frequency the localisation in R corresponds to a localisation in electron kinetic energy = 2 ( − 1) = 2 ( ⁄ − 1), where Rc is the cold resonance location, defined by ( ) = . The location of the maximum of the radiation function Frad in normalised momentum p/pth (where pth=mvth=(mTe) 1/2 ) depends on the harmonic number and on the electron temperature, as shown in Fig. 3, where the maximum is plotted together with the widths at half height. This figure illustrates a fundamental property of the ECE diagnostic: different harmonics probe different parts of the electron distribution function, both momentum and width decreasing with the electron temperature. This means that ECE at 2 nd and 3 rd harmonics can be used to constrain the electron distribution function in the region pth<p<2pth, in a more and more precise way for higher and higher temperatures.
Two examples of the previously defined isotropic model perturbation are shown in Fig. 4, together with the ranges seen by the 2 nd and 3 rd harmonics, for Te=7 keV.
It appears that the two harmonics can be differently affected by perturbations localised at different momenta, i.e., with a different p0 parameter. In general, the smaller the absolute value of the distribution function derivative in a given harmonic range, the larger the corresponding radiation temperature. Moreover, since the locations and widths of the momentum regions probed by the two harmonics vary with temperature, the same perturbation will affect the Te profile measured by ECE more or less significantly, depending on the temperature. An example is shown in Fig. 5, where the 2 nd harmonic ECE profiles are computed for typical JET parameters and two different values of the central temperature Te0: 3 and 10 keV. For the same perturbation parameters (f0 = 0.03, p0/pth = 1, /pth = 0.25), the radiation temperature is significantly affected at high temperature, but practically unaffected at low temperature, in agreement with experimental observations [1,2].
The strong sensitivity of ECE to very small perturbations of the electron distribution (a few percent in the example of Fig. 5) is due to the presence of the exponential term in Frad and to the fact that the absorption coefficient is an integral in momentum space of perpendicular derivative of the distribution function [10]. The same perturbation would have very little effect on the Thomson scattering measurement, which is simply proportional to the distribution function (see, e.g., Eqs. 5.8, 5.9 of [11]). This is illustrated in Fig. 6, where Frad is plotted as a function of the normalised electron kinetic energy, together with the equivalent quantity for Thomson scattering (scattered radiation spectrum). Clearly, the impact is completely different. This gives the main key to understand why the two measurements of the electron temperature can give different results if the distribution function is not Maxwellian.

Comparison with ECE and TS data
The distribution function model can now be used as an analysis tool of the extensive JET database described in  Ref. [9]. Examples of various subsets of the database are shown, corresponding to specific experimental scenarios.
In Figs. 7-9, ECE central temperatures measured at both the 2 nd and 3 rd harmonics by means of the JET Martin-Puplett interferometer [12] are shown vs the corresponding temperature measured by Thomson scattering via the so-called LIDAR system [13] (left panels). These temperatures are averaged over a region covering 10% of the minor radius around the centre. In the right panels, the central temperatures measured at the 3 rd harmonic are directly plotted vs the corresponding 2 nd harmonic ones. In Fig. 7 only, error bars are also shown. On both panels, the corresponding quantities computed using the above described model are also shown, choosing (by trial and error) a set of distribution function parameters that optimise the agreement with data, identical for both harmonics. The perturbation set of parameters used is displayed at the top of the figures. In all cases, a wall reflection coefficient of 0.55 has been assumed, according to previous evaluations for the JET machine [12]; however, the results are very weakly dependent on this parameter. The three figures refer to three different subsets of data, corresponding to different experimental scenarios [9]: - Fig. 7: DD baseline discharges with low gas, low Neon injection and pellets, characterised by small ELMs and partially detached divertor [14]. - Fig. 8: DT discharges in the baseline scenario. - Fig. 9: DT discharges in the hybrid scenario.  In all these figures, the 3 rd harmonic points are well below the corresponding TS measurements at low temperature (Te0 < 4 keV, approximately), because of the low optical thickness. Fig. 9. As in Fig. 8 In Figs. 7 and 8, it appears that the temperature measured by the 3 rd harmonic, in the optically thick range, is close to the LIDAR temperature, or slightly lower. Following the discussion of Figs. 3 and 4, this means that the distribution function is nearly unperturbed in the velocity range covered by the 3 rd harmonic (~ 2vth). On the other hand, the 2 nd harmonic temperature is significantly higher than the TS one at high Te0, meaning that some flattening of the distribution function takes place around vth. Indeed, the optimum fit is obtained for p0/pth = 1.1 -1.2. This situation is very similar to that of past JET experiments [2].  In Fig. 9 (DT hybrid scenario), the 3 rd harmonic tends to yield lower temperatures and the optimum fit is obtained for p0/pth = 1.4. There are also cases in which the 3 rd harmonic is strongly affected by the perturbation, whereas the 2 nd one is weakly affected. In general, phases with ICRH only have these characteristics. As an example, it is interesting to show how the model behaves for the simulation of an individual pulse (96850), characterised by low density and, in the high temperature phase, ICRH only. Figure 10 shows that in this case the 2 nd harmonic is nearly unperturbed, whereas the 3 rd harmonic has non-Maxwellian features that are well reproduced by a broad perturbation located at p0/pth = 2. The fact that ICRH tends to produce perturbations at somehow higher momenta is a general trend, that can be illustrated by considering all the DT discharges (belonging to both baseline and hybrid scenarios) and selecting data points with NBI only and with ICRH only. There is a significant number of them. Results for these two data sets are shown in Figs. 11 and 12, respectively. ICRH only phases of the DT discharges are clearly perturbed at higher momenta than NBI only ones, i.e., p0/pth = 1.7 with respect to 1.35. The intensity of the perturbation is also higher (f0 = 0.05 instead of 0.03). This is possibly related to the different types of ion tails that the two heating systems produce and/or to their different mechanisms and intensity of direct interaction with the electrons.
These examples have shown that the most sensitive parameter of the model is the perturbation location in momentum space, p0. This suggests a possible use of the model as an analysis tool, in order to detect trends in the database with respect to various distinctive quantities. This is illustrated in Fig. 13. Using the full database, the p0 value is determined by fitting the measured ECE central temperature (at both harmonics) for selected values of four quantities: central density, ratio of heating power and central density, Alfvén velocity normalised to electron thermal velocity and fast ion beta (this quantity is not a measured one, but is obtained from results of NBI and ICRH modelling codes that are available in the JET database). Slight adjustments of the other model parameters (perturbation width and intensity) are made in some cases. Figure 13 shows that a regular behaviour of p0 is obtained in all these cases; the trends observed can be used to guide the search for an interpretation of the experimentally observed non-Maxwellian features.

Conclusions
The discrepancy observed in JET plasmas [2,4,9] between ECE and TS measurements has been analysed with the hypothesis that it could be associated with a non-Maxwellian bulk electron distribution function. This hypothesis has been already formulated in the past [4,5], however, it remained to be quantitatively assessed on an extensive database. In order to perform such a quantitative comparison, a toy model of non-Maxwellian electron distribution function has been developed, allowing analytical computation of the Xmode ECE spectra and massive comparison with the JET central temperature database [9], with more than 12000 data points.
Using this model, various points have been clarified: an electron distribution function perturbation of a few percent, localised around 1-2 thermal velocities, is sufficient to explain the level of discrepancy observed; such a perturbation would be practically invisible to TS diagnostics, at least those of JET and other existing machines; however, specific TS systems can be conceived to this end [15]; even if such a perturbation is present at any temperature, it becomes visible and more and more significant at higher and higher temperatures, because the radiation function becomes narrower and shifts to lower velocities; at high temperature (Te0 > 4-5 keV), X-mode measurements at both 2 nd and 3 rd harmonic (possibly also higher harmonics and/or 1 st harmonic O-mode) are essential in order to properly constrain the distribution function in different velocity ranges.
An important open question is whether a perturbation of this kind could significantly affect measurements of the ECE profile in ITER. Of course, since the cause of the perturbation is at present still unknown, there is no reason to assume that the same mechanisms acting in JET (and TFTR) could also be present and significant in an ITER plasma. The model only allows quantifying whether a given perturbation would affect the ECE temperature profile in a sizeable way. This is what is shown by the example of Fig. 14. Two non-Maxwellian distributions are considered, with perturbations of the same intensity as those observed in JET plasmas, similar width and two different momentum localisations: p0/pth = 0.75 and p0/pth = 1. The impact on the 2 nd harmonic ECE profile is shown in the left panel (red and blue stars, respectively). It appears that at the high temperatures expected in ITER (25 keV in this example) the effect can be stronger and acting in two opposite ways, depending on the momentum localisation of the perturbation. Because of the large temperature variations on an ITER profile, both effects can be observed on the same ECE profile: note that in the 10-15 keV range, the blue ECE profile is higher than the Maxwellian, whereas it becomes lower beyond 15 keV. Therefore, it would not be just a matter of central temperature value, but a distortion of the full profile. The model described here should be regarded as an analysis tool only. It is not linked to a theoretical explanation about the origin of the electron distribution function perturbation. At present, we could just formulate hypotheses. Because of the nature of the JET high performance discharges analysed [9], which is analogous to those of TFTR [1], likely explanations involve the fast ion population, which is ubiquitous in these discharges (driven by NBI, ICRH and the fusion reactions themselves for DT pulses). In addition, ICRH also directly interacts with the electrons. Fast ions might interact with the electrons via two mechanisms: collisional relaxation of the fast ion tail on the electron distribution; interaction of fast ion driven MHD modes with the electrons (e.g., via Landau damping) About the first mechanism, it is generally assumed that the electron distribution function basically remains Maxwellian, because interaction with the fast ion tail takes place at low velocities (around thermal or subthermal), where collisions are very strong. However, we have demonstrated that a tiny perturbation (~a few percent) is sufficient to explain the observed effects. In this velocity range, asymptotic expansions of the collision operator cannot be used and, looking for small effects, even the usual linearization of the collision operator could be questionable. Therefore, a solution of the kinetic equation with the full integro-differential collision operator is likely to be needed, which is not generally available in the literature. In Ref. [6], this type of problem is solved for inertial fusion applications and the resulting perturbation has a clear bipolar structure, in the thermal velocity range. Work is now ongoing in order to solve the problem numerically in the JET parameter range and for various input fast ion distributions (NBI, ICRH driven tails or alpha particle distributions can be computed by means of Monte-Carlo or Fokker-Planck codes).
About the second mechanism, direct observation by probes in the magnetosheath have revealed the presence of bipolar distortions of the electron distribution function, around the electron thermal velocity [7]. Gyrokinetic simulations have provided a convincing interpretation of such observations in terms of Landau damping of kinetic Alfvén waves [8], which are also known to be present in tokamak plasmas. Nevertheless, magnetosheath plasmas have very different characteristics with respect to tokamak plasmas (for instance, they are nearly collisionless), therefore, these results cannot be easily generalised to the case of interest here. However, they have inspired analogous gyrokinetic calculations that are now in progress, for parameters close to those of the JET experiments. Note that, in principle, other MHD modes (not necessarily directly excited by fast ions) could also interact with the electrons and cause distribution function distortions.
In conclusion, as for the diagnosis of the electron temperature profile in high temperature plasmas, TS appears rather insensitive to small perturbations of the electron distribution function, therefore, in this respect, it is expected to provide a reliable measurement of the electron temperature. Conversely, ECE can definitely be affected by tiny perturbations of the electron distribution (a few percent) localised in the range 1-2 vth, in a different way at different harmonics and in different temperature ranges. Therefore, the temperature measurement can be considered as less reliable than that of TS. However, this high sensitivity of ECE can be exploited to constrain the electron distribution function in order to extract information on its detailed shape and explore fundamental physics effects, such as, e.g., those related to fast ion physics.