J Vicente
email: jvicente@ipfn.tecnico.ulisboa.pt

T Ribeiro

F Da Silva

S Heuraux

G D Conway

C Silva

2D full-wave simulations of conventional reflectometry using 3D gyro-fluid plasma turbulence

Keywords: reflectometry, turbulence, full-wave, gyro-fluid, simulations measurements

Asynthetic reflectometer based on the 2D finite-difference time-domain full-wave code REFMUL is applied to data from 3D numerical simulations of edge and scrape-off layer plasma turbulence obtained with the GEMR code. Full-wave simulations are performed using the conventional reflectometry set-up with O-mode waves, fixed frequency probing and an equivalent I/Q detection scheme. Results show a significant spectral broadening of the synthetic reflectometry complex amplitude A(t)e ij (t) with increasing probing frequency. The reflectometry response displays 2p phase jumps which are due to the self-consistent evolution of turbulent density structures. The range ∼3%-6% of moderate turbulence amplitude is studied here, in the transition from the linear to nonlinear regimes of conventional reflectometry. While a phase jump removal algorithm is applied, spectral broadening of the phase with increasing probing frequency is nevertheless observed. Linear scaling of phase fluctuations with d /nn ee is also retrieved. REFMUL simulations with turbulence data rendered on both the GEMR field aligned coordinates (drift planes, neglecting circular magnetic flux surfaces) and on polar geometry (poloidal planes, taking into account plasma curvature) are carried out revealing similar trends on both coordinate systems. However, phase fluctuations obtained from poloidal planes display higher root mean square values, compared to drift planes. This could be expected from an increased sensitivity of the reflectometer to higher wavenumbers, due to plasma curvature effects.

Introduction

The performance of magnetized plasmas in fusion devices is widely determined by turbulence and magneto-hydrodynamic instabilities. To understand these phenomena, diagnostics with high temporal and spatial resolution must be employed to allow experimental characterization of the relevant plasma properties. Microwave reflectometry techniques have been successfully applied to this end, measuring for instance electron density fluctuations (see [1] for an overview). Nevertheless, without a strong modelling effort, turbulence turbulence. Recently, complete chains from measured turbulence properties through a full-wave code simulating reflectometry in realistic gyro-kinetic simulations have been implemented [START_REF] Da Silva | Proc. 9th Int. Reflectometry Workshop -IRW9[END_REF][START_REF] Hacquin | [END_REF][8][9]. Modelling of edge and scrape-off layer (SOL) plasma dynamics with a fluid model equipped with other synthetic diagnostics has also been reported [10]. In this work, the two-dimensional finite-difference time-domain (FDTD) full-wave code REFMUL [11] is used to implement synthetic conventional reflectometry diagnostics on turbulence simulations carried out with the gyro-fluid code GEMR. In particular, the electron density dynamics are diagnosed with fixed frequency O-mode reflectometry in the plasma periphery, including the plasma edge and SOL. In this region, typical experimental background density profiles and density fluctuation levels / dn n e e are as displayed in figure 1 [12,13]. This figure illustrates important features of the underlying physics and their impact on the reflectometry technique, as will be shown later. Note, for instance, the steeper density profile region in the vicinity of the Last Closed Flux Surface LCFS, and the drastic change in the fluctuation levels from core (<1%) to SOL (>20%) where intermittent high density perturbations are dominant. In this work, the edge to SOL transition (shaded region in figure 1) is the region of interest where fluctuation levels increase significantly and more challenging reflectometry measurements will be investigated. In addition, the impact of plasma curvature on reflectometry measurements will also be assessed together with the validity of the usual assumption of taking data rendered in the socalled drift plane equivalently to the poloidal plane. The paper is organized as follows: in section 2 an overview of the technique and theory of conventional reflectometry for fluctuation measurements is provided. The turbulence code GEMR and the full-wave code REFMUL are briefly described in sections 3 and 4, respectively. In section 5, the main issues and adopted strategies for coupling the two codes are explained. Details about the set-up used for both GEMR and REFMUL simulations are given in section 6. In section 7 the procedure to obtain synthetic reflectometry signals from REFMUL in the frozen turbulence approximation is given. In section 8 a preliminary characterization of GEMR output turbulence is made without analysing synthetic reflectometry data. In section 9 the results of synthetic reflectometry are presented and compared to turbulence characteristics and theoretical predictions given previously. Finally, a summary and brief discussion of the results is done in section 10.

Conventional microwave reflectometry

Conventional reflectometry is based on the phase variation of a microwave that is launched into the plasma, with incidence along the plasma surface normal, and reflected at a cut-off layer where the refractive index N becomes zero. As O-mode polarized waves with fixed vacuum probing frequency f o propagate in the plasma, its wavenumber decreases until the cut-off condition is met, and reflection occurs, at a critical density n c given by ⎛ where e and m e are the electron charge and mass, and e 0 is the permeability and permittivity of free space, respectively. Fixed frequency waves can thus probe different regions of the plasma according to the plasma radial density profile.

Early developments of wave scattering and microwave reflectometry theory [14][15][16][17][18] supported the development of reflectometers aiming to measure density profiles and density fluctuations [19][20][21][22]. The phase of O-mode probing waves, propagating in the plasma, records the reflection at the cut-off but also the fluctuations due to plasma density perturbations. This includes Bragg backscattering processes, which are plasma-wave interactions of the incident wave with density perturbations satisfying the Bragg rule ( )

   = - k k r 2 f i B
where   k f is a given wavenumber component of a density perturbation and  k i is the local wavenumber of the incident wave at the radial location r B along its path in the plasma. Density fluctuations displaying only radial wave vector   k r are considered in light of one dimensional models, for which the Bragg resonance rule can be satisfied in the range

< < k k k 2 2 , A r o
where k o is the vacuum probing wavenumber and k A is the Airy wavenumber 11 3 with L n being the density gradient scale length [23]. The minimum plasma fluctuation wavenumber allowing Bragg backscattering, = k k 2 , r A resonates in the vicinity of the cut-off. Under the Born approximation, for small amplitude perturbations and a slab background plasma, a linear relation between the O-mode phase and the local density perturbation level / dn n c has been shown [23,24]:

( ) / = - k k L 0.63 , A o n 2
⎛ ⎝ ⎜ ⎞ ⎠ ⎟ ( ) dj p d = k L k n n . 2 o n f c
A common 1D result for fusion plasmas was the localised phase response from the reflecting layer and fluctuations in its vicinity [25]. To obtain a realistic description of experiments, two dimensional effects (e.g. refraction and finite beam divergence) should be considered even if neglecting any density fluctuations. To analyse the phase response in 2D geometry, various models have been proposed including the physical optics approach [26,27] and full-wave solutions [28][29][30]. It was shown that the reflectometer was sensitive to poloidal fluctuation wavenumbers (poloidal and perpendicular directions assumed identical) k up to a limit determined by the receiver location and incident beam width (or spot size). A sensitivity criterion < k w 2, where w is the 1/e radius of the incident beam intensity, was deduced from rigorous analysis based on slab geometry [31]. Later it was found that the plasma poloidal curvature and the curvature of the probing wave-front could extend the microwave reflectometry response to much higher poloidal wavenumbers than what was prescribed by slab plasma models [32][33][34]. A phase screen model was used to derive the strong curvature regime criterion

( ) /  r k w 2 o 2 1
, where r is an effective radius of curvature defined by [32]:

( ) ( ) / r r r r r = + 2 2 3
c w c w using the radius of curvature of the cut-off layer r c and of the incident wave front at the cut-off layer r . w Under strong curvatures, the sensitivity to fluctuations with a given wavenumber k showed to be as strong as to much smaller wavenumber ¢ k without considering curvature effects:

⎛ ⎝ ⎜ ⎞ ⎠ ⎟ ( ) r ¢ = ^k k w k 2 . 4 o 2
Note that both the beam width and the wave front curvature depend on the propagation in the plasma (e.g. [34]). On the other hand, the interpretation of correlation reflectometry for measurements of fluctuation levels and turbulence correlation lengths proved to be challenging. Detailed knowledge of the instrument transfer function was also shown to be required, for instance depending on which reflectometry signal would be used, such as the complex amplitude, the phase or the homodyne signal [35][36][37].

Forward scattering or small-angle scattering off longscale fluctuations ( » k 0 f ) may also contribute to reflectometry fluctuations, either considering 1D [17] or 2D models [38][39][40]. While these effects can become important, in particular under high turbulence conditions or in long plasma paths with stronger spreading of the probing beam [41], at low turbulence levels and small to mid-sized devices, Bragg backscattering in the vicinity of the cut-off is assumed to dominate. In any case, qualitative measurements are obtained mostly in the frame of the so-called linear (scattering) regime of reflectometry, which is in practice limited to small phase fluctuations and turbulence levels [25,38,42,43]. A transition to the nonlinear regime of reflectometry has been explained via small angle multi-scattering mechanisms [38]. It was concluded that in the nonlinear regime the reflectometry spectrum is not localized to the cut-off vicinity conveying information about a wide region of the plasma, despite a quick coherence decay remaining proportional to the local level of density fluctuations at the cut-off. A criterion for the transition to the nonlinear regime was derived using a 1D modelling with linear density profile and statistically homogeneous turbulence with arbitrary wavenumber spectrum:

⎛ ⎝ ⎜ ⎞ ⎠ ⎟ ( ) x x d w =  n n l x c x l 1, ln . 5 c c c c c 2 2 2
This criterion, where l c is the turbulence correlation length, w is the angular wave frequency, x c is the plasma length to cut-off and c is the speed of light in vacuum, showed to be valid when 2D effects were also considered in 2D fullwave simulations [44,45].

With increasing density fluctuation levels, phase jumps start to occur [25,29,46]. If predominant, phase jumps lead to distorted reflectometer spectra, approaching a f -2 decay over most of the dynamic range of the diagnostics [20,47,48]. For this reason, removal of phase jumps has been carried out in experimental data mostly on an empirical basis [49]. A convolution method was suggested [50] and allowed to retrieve both the spectra and correlation lengths of turbulence in good agreement with probe measurements [51]. Other authors have employed either similar algorithms [12,52,53], back-projection techniques [25], or simply bandpass filtering to cope with phase jump and phase runaway phenomena. If only a small fraction of the phase data corresponds to phase jumps this may be indicative of a transition in the reflectometer response regime. Independently of the method employed, phase jump removal may be a reasonable way to extend the range (only up to a certain limit) over which phase measurements might be used in linearly determining plasma turbulence levels.

A fully realistic description of the phase jump phenomena may require employing 3D modelling in order to properly describe the behaviour of electric fields in 3D cavities. REFMUL3, a three-dimensional full-wave code with parallel computing, was just recently developed and first benchmarked against the 2D code REFMUL [54]. It was shown a stronger impact of 3D effects in the reflectometry amplitude rather than in the phase for slab plasmas. Although this first result was obtained without plasma turbulence, it strengthens our confidence in 2D simulations with respect to describing reflectometry phase robustly, while 3D simulations are still computationally expensive and a scarce resource.

Turbulence code: GEMR

The turbulence code GEMR consists of a three-dimensional (3D) electromagnetic gyrofluid model with global geometry suitable for edge plasma conditions. The code solves the first six moments of the gyrokinetic equation, computing the evolution of density fluctuations and the density profile gradient. To yield a closed system, the gyrokinetic polarisation equation [55] and the Ampère's law are also considered. The gyrofluid formulation guarantees that finite Larmor radius effects (FLR) are included in the model. In fact, whenever ion and electron temperatures are T i ~Te , and hence the ion gyroradius (r i ) and the size of the typical vortical space scales (r s ) are of the same order r r ~, i s the dynamics on the scales r k i 1 become important. This is relevant for this work since the spatial resolution demanded by the reflectometer to produce accurate phase signals goes beyond this range (see section 4).

The GEMR geometry is global in the sense that both poloidal and radial dependencies of the geometrical quantities are considered, without using the flux tube approximation. The GEMR coordinates are field aligned, hence non-orthogonal, which is justified by the computational efficiency gain allowed by the strong spatial anisotropy between the direction along the magnetic field and the plane that is perpendicular. A generalised poloidal angle is set as the field aligned coordinate ϑ, while ξ follows the toroidal angle to give the field line label and ψ is the (radial) flux label. The simulations presented in this paper use a circular equilibrium described according to the unit Jacobian definition

( ) y p = R r 2 , 6.1 2 0 2 ( ) ( ( ) ( )) ( ) x p q f a = -- -q r r 2 , 6 . 2 k k 1 ( ) ( ) J p q = - 2 , 6.3 1
where ( ) q f r, , are the polar coordinates, R 0 is the tokamak major radius and ( ) q r is the magnetic field pitch. The quantity ( ) a r k represents a toroidal shift due to the shearing of the magnetic field, which is applied locally along the magnetic field lines to ensure a correct representation of both toroidal and slab mode structure. The details of such technique called shifted metric procedure can be found elsewhere [56].

An example of a representative GEMR turbulence snapshot for a medium-sized tokamak, like ASDEX Upgrade (AUG), is shown in figure 2. The very low ratio between the radial and poloidal domain sizes can be noted in this case. Also displayed in the figure is a zoom-in view of a typical region of interest on the equatorial low-field side (LFS).

Because the background profiles evolve in time, a dynamical Shafranov shift and a correction to the magnetic field pitch are calculated to adjust the field alignment self consistently in the model [57]. Global consistency in the angles is also ensured [58], even in toroidally truncated domains, such as in the cases presented here, which simulate one fourth of the whole torus. This is a major point in the closed field lines region, where it is insured that parallel space scales follow the so called field line connection constraint

 k qR 1. 0
In the SOL the field lines cease to wrap around the flux surface to intersect material plates, allowing the existence of modes with both  = k 0 and ¹ k 0 (finite-sized eddies in the perpendicular plane which are perfectly aligned to the magnetic field also called convective cell modes,). The fact that convective cell modes are allowed outside the LCFS, and not inside due to the field line connection constraint, changes the character of turbulence from drift wave to an interchange turbulence [59]. In gradient driven turbulence of magnetised plasmas, the balance between linear and nonlinear processes determine the evolution of the system. Nonlinear processes include an inverse energy cascade due to vorticity nonlinearity [60] and a direct cascade due to E×B advection of pressure. Among linear processes are the interchange forcing and the so called adiabatic response where parallel electron dynamics act on all available degrees of freedom, keeping the electrostatic potential ( f) coupled to the electron pressure ( pe ) through parallel forces mediated by parallel currents (  J ) [61,62]. Since convective cells do not experience parallel dynamics, in the SOL plasma f and pe are coupled through the interchange mechanism which becomes the dominant turbulence character. Closed and open field line regions thus experience substantially different behaviour of heat and particle fluxes down the gradients, as previously demonstrated with GEMR [59].

Full-wave code: REFMUL

REFMUL is a full-wave propagation code that solves the Maxwell equations in two dimensions for ordinary wave polarization (O-mode) using the FDTD technique. The plasma electron density ( ) n r t , e sufficiently accounts for the plasma effects since O-mode propagation fulfils the cold plasma approximation. Both electric and magnetic fields of the electromagnetic wave are coupled to ( ) n r t , e in the Yee scheme through the density current J [63]. The complete set of FDTD equations is the following

( ) ( ) ( ) m e m ¶ = ¶ - ¶ - E B B J , 7. 1 t z i j n x y y x z i j n 0 0 , 0 ,
) ( ) ( ) / / / / ¶ = - ¶ + + + + B E , 7 . 2 t x i j n y z i j n , 1 2 1 2 , 1 2 1 2 ( ) ( ) ( ) / / / / ¶ = - ¶ + + + + B E , 7 . 3 t y i j n x z i j n 1 2, 1 2 1 2, 1 2
( ) ( ) ( ) / / / ¶ = + + J e m nE , 7 . 4 t z i j n e e z i j n , 1 2 2 , 1 2
where m 0 is the permeability of free space and the indices ( ) i j , refer to the spatial coordinates. The term n e n provides the plasma model (at the time instance = ´D t n t) which in this work is obtained from results of GEMR simulations. A detailed description of GEMR data processing, required for REFMUL integration, will be given later in section 5. The REFMUL code runs on a Cartesian grid of dimensions L x ×L y where the x-and y-directions correspond to the radial and poloidal directions in a tokamak. Part of the simulation box is devoted to the signal emission setup (antenna) and a vacuum distance to the plasma. A unidirectional transparent source is used for injection of the signal in a metallic structure, as a waveguide or an antenna, allowing separating the emitted probing wave from any returning waves [11]. A Perfectly Matched Layer (PML) approach is used for boundary conditions [64]. In this work, the response of O-mode reflectometry in conventional set-up was studied using a monostatic design, i.e. one antenna for both emission and reception.

Separate REFMUL runs with different fixed frequency values have been performed, in particular to allow probing of regions inside and outside the LCFS. To be as close as possible to the experiment, the same reflectometer parameters are extracted using similar data processing tools. Moreover, the synthetic reflectometry phase j(t) and amplitude A(t) signals are obtained carrying out identical simulation runs with sin and cos excitation signals allowing in-phase and quadrature detection, as will be explained in section 7. To adequately resolve the phase j(t), the spatial discretization in FDTD simulations must be a small fraction of the probing wavelength. In this work, simulations are run with

/ d d l = = = ´- x y 20 3.75 10 m, 40GHz 4
where l 40 GHz is the vacuum wavelength of a wave with frequency f=40 GHz which could still be suitably employed. The time discretization is also a small fraction of the wave period in order to comply with the Courant-Friedrichs-Lewy condition to ensure stability. The temporal resolution is thus required to be d =

´t 6.25 10 s. 13 This makes reflectometry simulations computationally very demanding, which is the main reason why 2D codes have been used rather than 3D.

Coupling of GEMR and REFMUL

Different codes use distinct models, approximations and mathematical descriptions. A direct coupling between them is often challenging. In particular, to integrate GEMR and REFMUL special attention must be paid to the coordinate systems employed in each code and the required spatial and temporal resolutions.

Poloidal plane and drift plane

The main approach in this work is to use instances (2D maps) of plasma electron density provided by GEMR as the required plasma model input to REFMUL simulations. The electron density on the drift plane ( ) y x , perpendicular to the magnetic field as defined by the GEMR coordinates can be used without further coordinate transformations. This approximates the drift plane with the poloidal plane and neglects the circular cross section of the magnetic flux surfaces. Previous methods to couple numerical turbulence simulations and reflectometer simulations have employed the drift plane approximation (e.g. [2]). However, experimental reflectometry probes the plasma electron density on the poloidal plane where plasma curvature effects may play an important role, as highlighted in section 2. For that reason, post processing coordinates transformations to invert equations equations (6.1)-( 6.3) and express the electron density in the usual cylindrical coordinates has to be done. This involves transforming into Fourier space in the GEMR toroidal angle coordinate (ξ) to apply a phase shift that undoes the shifted metric procedure [56]. This must be followed by interpolation from a typical lower resolution parallel grid to a high resolution one, after which a phase factor of q (magnetic field pitch) is applied to transform back to an unaligned representation-that is, to move from a parallel coordinate to a poloidal coordinate ( ) q . To finalise, an inverse Fourier transform is applied by summing the transformed toroidal Fourier modes, which yields a single poloidal plane where the data n e is then measured. In summary, we employ both the poloidal plane and the drift plane rendering of the n e data on REFMUL simulations in order to compare the results obtained with both methods, in particular to assess any effects of the plasma curvature on the reflectometer sensitivity.

Interpolation

The time discretization usually employed in GEMR is much larger than in REFMUL. A GEMR time resolution of d ¢ =

´t 9.3 10 s, 7 allows to consider the plasma frozen in the time frame of the probing signal.

As previously mentioned, electron density data obtained from GEMR can be provided either on a drift plane cut or on a poloidal plane cut of the plasma. The former is defined on a Cartesian grid, as expected for REFMUL, and interpolation is done merely to ensure the required spatial resolution (number of points per wavelength) is met. The latter is defined on polar geometry ( ) q R, and the points of the two meshes do not obviously coincide. Furthermore, the points of the polar grid are not equidistantly distributed. This problem is treated with the usual tools for irregular meshes, such as Delaunay interpolation, in order to obtain the required spatial resolution. In practice, the interpolation is not performed on the whole poloidal cut, but only in a smaller region, adequate for the REFMUL simulations. The entire poloidal plane cut domain of an AUG turbulence snapshot and a typical region of interest on the LFS were shown previously in figure 2.

Boundary conditioning

One must be aware that the GEMR density profiles do not start at zero density. The numerical SOL in these particular GEMR runs has a minimum/outer value » ´- n 9 10 m e 18 3

which is relatively high. Employing probing waves on such density profiles with frequencies near the value corresponding to that critical density ( º f o 27 GHz, according to equation (1)) would result in sudden phase jumps as the waves would propagate from vacuum suddenly into the plasma medium. To overcome this, the numerical plasma is extended to allow a smooth growth from ~0 density to the first GEMR density point. A modified hyperbolic tangent [65] is used to shape this numerical SOL extension which has a length of 80 grid points. Similarly, to allow transitioning from the plasma medium to the perfectly matched layer used for boundary conditions, the plasma is extended towards the other edges of the simulation grid by 100-150 grid points to accommodate adaptation layers. In figure 3 a radial cut (profile) taken at the middle poloidal position on the LFS is shown for a Cartesian density matrix prepared to be used as input for REFMUL simulations.

Simulation Setup

The GEMR simulations assumed a simplified magnetic equilibrium with circular plasma cross-section and local plasma parameters representative of L-mode discharges at the AUG tokamak (e.g. see [66] for a comprehensive range of ohmic L-mode discharges used for intrinsic rotation studies). Namely, the following parameters were used: major radius R=1.65 m, minor radius a=0.5 m, magnetic field B t =2.6 T. The plasma characteristics are ion and electron temperatures T i =T e =60 eV, ion and electron densities

= = ´- n n 1.2 10 m i e 19
3 respectively and temperature and density gradient scale lengths L T =L n /2=3.5 cm, nuclear mass of deuterium M D =3670m e , effective nuclear charge number Z eff =2, safety factor q=4.6, and magnetic shear ˆ= s 1.13. In this work, we focus on the role of edge and SOL turbulence on the synthetic reflectometry measurements, which requires GEMR to simulate a somewhat broad radial plasma region that includes regions with closed and open flux surfaces.

The GEMR simulations were run allowing the evolution of turbulence using (ψ, ϑ, ξ)=(512×1024×16) grid points and recording electron density matrices (drift plane cuts and poloidal plane cuts at the mid-plane LFS) at 1915 sequential time iteration points. Subsequently, each GEMR matrix, was processed as explained in section 4. Both drift plane and poloidal plane density matrices could then be used as input to REFMUL reflectometry simulations. Independently of the type of plane cut, the plasma itself was recorded on a similar region with dimensions l x ∼ 6 cm and l y ∼ 17 cm. The whole REFMUL simulation box was larger, due to the boundary extensions previously described, as well as additional room for the emission setup (antenna) and a vacuum distance (d vac =413 grid points) to the plasma. In this work a monostatic design was employed (i.e. one antenna for both beam emission and reception) together with a unidirectional transparent source, allowing separation of the emitted probing wave from any returning waves. A 2D H-plane horn antenna with a half power beam width 6.5 cm at the plasma entry was chosen for directivity, with on-axis line of sight at the LFS mid-plane. The same antenna was employed in another similar set-up, where a lens focused the beam ahead of the cut-off positions, resulting in the half power beam width 3.2 cm.

O-mode probing waves with fixed frequency f o ={28, 29, 30, 31, 32} GHz, corresponding to cut-off densities in the range [0.97-1.27] ×10 19 m -3 , have been used. These density layers correspond to radial locations across the open (SOL) and closed (edge) field line regions, away from the boundaries where GEMR source and sinks act. The chosen set of frequencies also lie on fairly constant background gradients to ensure some independence of the reflectometer response to this parameter. Frequency values lower than 28 GHz would also yield smaller gradients in the outer regions, possibly rendering poor radial probing localization.

Synthetic reflectometry probing

Simultaneous access to the simulated GEMR density profiles and density fluctuations, allows to easily strip the latter from the former by averaging the data over the angles (zonal average). Hence, for each of the GEMR renditions of density (drift plane and poloidal plane), two sets of sequential n e GEMR data matrices were obtained and used as input to REFMUL. One set comprised the zonal averaged plasma density profiles ( ) n x y , (turbulent plasma). Probing of the two sets, separately, can be combined in order to implement a synthetic in-phase and quadrature (I/Q) detection that allows to compute the returned amplitude ( ) A t and phase ( ) j t reflectometry signals. This is performed as follows: two REFMUL runs are carried out per matrix, one with a cosine (I) excitation signal and the other with sine excitation (Q). The in-phase and quadrature emitted signals are, respectively: Finally, yielding the phase variations and amplitude as:

( ) ( ) ( ) p = s t
[ ] ( [ ] ) ( [ ] ) ( ) dj dj dj = - n n n tan sin cos , 13.1 1 [ ] ( [ ]) ( [ ]) ( ) dj dj = + A n n n sin cos . 13.2 2 2
Note that [ ] j n p can be obtained using the same procedure as above by similarly sampling the emitted signals (equation ( 8)) and combining them with equation (11) time series. In this work, and for the results presented in the following sections, the REFMUL runs were carried out with 12 0000 time iterations and the reflectometry phase and amplitude signals were obtained by sampling at time t k corresponding to n = 90 000, which was observed to be well into the stationary regime of the reflectometer response.

Characteristics of GEMR turbulence

Plasma density fluctuations contain a wealth of information about the underlying turbulence mechanisms [67]. Their perpendicular wavenumber spectra can be usually measured with Doppler reflectometry diagnostics, e.g. see [68]. Efforts to understand both the experimental and the theoretical results have been made, in particular employing synthetic diagnostics for validation studies [69]. Conventional reflectometry may also allow to retrieve the radial wavenumber spectra [70], however it is usually relevant to simply obtain the turbulence level of the electron density fluctuations.

The typical scale length of turbulence dynamics is given by the correlation length, defined as the distance by which the correlation between fluctuating quantities drops from its maximum value by a factor of 1/e. Cross-correlation functions (CCF) can be calculated between time series of electron density evolution registered at chosen grid points. Selecting a given poloidal position and computing the CCF for varying radial separations, allows estimating the radial correlation length l c directly from the radial separation value where CCF ( t D = r, 0)=1/e, where t is the time lag [33]. In figure 4, CCFs calculated this way using GEMR data, are shown for a set of radial locations in the mid-plane.

The reference radial locations employed in the CCFs were chosen to correspond to the cut-off layer positions, in the base plasma, of the set of probing frequencies to be employed in REFMUL (see figure 3). The radial correlation lengths l c obtained in the case of drift plane data, shown in the figure, were in the range of 0.31-0.44 cm. Poloidal plane cut data analysis resulted in similar values.

To further assess the local effects of turbulence on GEMR data, we have chosen to register the time sequences of the average density n e at those different cut-off locations, where the averaging was defined in a given spot size region. The spot sizes were assumed rectangular with radial length w R and poloidal length (in poloidal plane cases) or perpendicular length (in drift plane cases) w .

Z The radial length was The computed time traces of average density within each spot size as well as the corresponding spectra are shown in figure 6, for the case of drift plane data. The spectra obtained from the fluctuating components at all positions are rather similar, both in drift planes and poloidal planes (the latter are not shown). The typical spectrum is relatively flat at low frequencies and then decays with a given index. In particular, there is no significant change in the spectral index with f o but only a slight increase of power with increasing probing frequency, which is observed across the whole frequency range.

These results indicate a small radial dependence of turbulence effects at the cut-off, at least measured with such spatial averaging. The radial coverage included locations inside and outside the LCFS, where drift waves and interchange instabilities dominate, respectively. But, in fact, no clear changes are expected in the spectra of density time traces registered in the two regions [71].

Despite this result being expected, to investigate the sensitivity of the spatial averaging method to variations in the radial and poloidal/perpendicular width of the spot size, scans on these parameters were also performed. A scan of the perpendicular width of the spot size is shown in figure 7, while a scan of the radial width of the spot size is shown in figure 8. The scans used the case f o =28 GHz as reference, with = w 0.94 cm R and = w 12.4 cm Z which was then varied by employing a multiplication factor on one of the spot size dimensions while the other remained constant. The total power across the whole frequency range, but in particular in the low frequency range, is observed to increase with decreasing w Z in figure 7, but the spectral index remains unchanged across the whole scan. On the other hand, it is observed in figure 8 that variations in the radial size w R have a small impact on the low frequency range. However, reducing the radial size results in a shallower spectral decay, i.e. smaller spectral index.

Radial profiles of turbulence level can also be computed directly from GEMR data. The antenna on-axis line of sight was considered to calculate the radial profiles displayed in figures 9 and 10 corresponding to the drift and poloidal plane cuts, respectively. In both figures are shown the base density profiles together with the standard deviation of the density fluctuation level, registered at each radial position for the full time evolution of the simulation run. The two turbulence level profiles, from the two plane cuts, are quite narrow, nonmonotonic and decreasing outwards, across the LCFS. This is not in agreement with what is observed in experiments, where in the edge and near-SOL regions there is an increase, or perhaps flattening, of the turbulence level towards the plasma periphery (remember figure 1) and typically reaching much higher levels in the SOL.

The reason for this discrepancy may result from GEMR iterating self-sustained base profiles, which are solved for ( ) n ln e rather than n .

e In practice, an average value of n e is assumed to be constant while its gradient is evolving. However, this assumption cannot be realistic for edge and SOL regions. Hence, while GEMR generated turbulence may be representative in the whole region of interest, the base profiles may not, which is also manifest on the / dn n e e profiles. As future work, one may consider taking the fluctuation density data obtained here together with another base plasma data (custom building steeper base plasmas, or using experimental profile data) to obtain / dn n e e profiles closer to what is observed experimentally in the edge and SOL. In the scope of this work, it is worth pointing out that only mild differences between the two plane cuts are observed, and / dn n e e values at the cut-off locations of the considered probing frequencies are in the range of ∼[3-6] %. Visual inspection of this data immediately reveals distinct intermittent characteristics in the two signals. The complex amplitude signals

( ) ( ) ( ) ( ) = + j A t e I t iQ t i t
were also obtained and their spectra analysed. This is shown in figure 12, where the full set of probing frequencies, i.e. f o ={28, 29, 30, 31, 32} GHz, are represented for the drift plane cases. The frequency spectra have been renormalized to / f f o in order to account for the increasing Bragg backscattering efficiency with increasing frequency. In general, the spectra are fairly symmetric. The typical spectrum shows the DC component at zero frequency, accounting for the reflected signal, and the spectral broadening around the central component, accounting for density fluctuations at the cut-off and Bragg resonant backscattering along the plasma path [24,70]. Data in both drift plane cases and poloidal plane cases (not shown here), display a trend of spectral broadening with increasing frequency. However, the DC reflected power appears to be lower in poloidal plane cases. The same trends were obtained with the focused antenna set-up, yet displaying significantly higher levels of power across the whole frequency range (see figure 13 for an example).

The spectral broadening effect appears to be independent of the local turbulence level since the latter does not vary monotonically in either poloidal or drift plane cases (see figures 9 and 10, respectively). Increased Bragg backscattering might still be determining the spectral broadening. Either due to geometrical effects, such as beam spreading, allowing poloidal components contributions [72] or due to longer wave paths into the turbulent plasmas enabling contributions from locations farther away from the cut-off.

The total integrated power of each spectrum, normalized to its DC power component, is shown in figure 14 against both the turbulence level and the plasma path length, computed at the corresponding cut-off layer positions in the base plasma. The root mean square of turbulent density fluctuations can be estimated through the squared spectral power [73]. However, in both planes, a strong linear relation with the plasma path length is also suggested by figure 14. The two dependencies are likely entangled on the total integrated power response of reflectometry. Note that the power integration strongly depends on the low frequency range components (DC included), which are usually lower with increasing probing frequency. Without the DC normalization of each spectrum, data would be more scattered without easily distinguishable trends.

The time evolution of the reflectometry phase can also be obtained, as previously explained in section 6. An example of base plasma phase ( ) j t p and fluctuating phase ( ) dj t A t e i t for a set of probing frequencies f o using the drift plane cut and the poloidal plane cut, with and without a focused antenna. Spectra have been normalized to the DC components.

computed for each probing frequency employed on drift plane simulations is shown in figure 15.

Rapid phase variations are observed in ( ) dj t while ( ) j t p only varies slowly, as expected. The smooth variations of the zonal density profile (base plasma) are not responsible for the phase jump effects, as confirmed by the evolution of ( ) dj t in the right panel of figure 16. Turbulent density fluctuations cause phase jumps into another position in the phase plane. The phase jumps are not related to generic phase unwrapping issues (for instance due to insufficient sampling rates) but are a common occurrence in experimental data, as discussed previously in section 2 (e.g. see [74]). Simulations have been designed previously to understand this phenomenology, but have always been done by imposing simplified models [46,75]. In the present work, the self-consistent evolution of the turbulent plasma forms the required conditions for phase jumps, being in this sense much nearer to reality and the experimental observations. To further perform spectral and statistical analysis on the phase signals, the phase jumps have been removed using an algorithm based on the one suggested by Ejiri [50]. An illustrative example is shown in figure 16 where the synthetic phase components obtained from REFMUL are displayed on the left and the corrected fluctuating phase signal is displayed on the right.

The resulting phase fluctuating signal components, obtained after phase jump corrections, are shown in figures 17 and 18, for the drift plane and poloidal plane respectively, as well as the corresponding power spectra. In the low frequency range, the shape of the spectra is relatively flat, in both drift and poloidal planes, and there is a strong roll-over of the spectra at a given knee position, followed by a region of constant spectral index, in agreement with the shape of the spectra previously obtained from the GEMR turbulence directly (see figure 6). A spectral broadening effect, observed previously in the complex amplitude signals, is also present in the phase signals. The steep spectra become increasingly shallower with increasing probing frequency, reaching higher spectral index values than what was obtained directly from the 2D density averaging procedure. At some probing frequencies ( f o =30 GHz and f o =31 GHz) there is clearly a larger power response at high frequencies in the poloidal plane simulations than in drift plane cases. This could be an effect of the plasma curvature, accounted for by the base plasma data in the poloidal plane cuts, which has been reported to extend the microwave response to high perpendicular wavenumber fluctuations [32].

However, in both planes, the differences in the spectra obtained at different probing frequencies (see figures 12, 17 and 18) are much stronger than those observed with the simple local 2D density sampling procedure (see figure 6). The latter, of course, did not account for any Bragg resonant effects. Differences in the dimensions of the spot sizes, alone, could not result in the large spectral modifications that are observed with synthetic reflectometry. This highlights the sensitivity of the reflectometer response to the beam broadening and/or increased Bragg backscattering efficiency away [76]. Using the radial correlation lengths computed previously in section 8, we obtain at the corresponding cut-off locations k f =3.2 cm -1 for f o =28 GHz, and k f =2.3 cm -1 for f o =32 GHz.

The theoretical estimates for the lowest and highest probing frequency draw a theoretical envelope for phase fluctuations due to Bragg backscattering. The poloidal pane data is mostly within the predictions given by equation (2). Note that only considering higher effective turbulence wavenumbers k f would bring the theoretical curves below the synthetic data, which would be unreasonable, at least for the lowest probing frequency. This effect might be due to the phase jump correction carried out for the poloidal plane cases, given the existence of at least one phase jump in most phase time traces. On the other hand, the phase fluctuations in the drift plane cuts are over-estimated by the theoretical model. In the drift plane cases only the phase at the highest probing frequency displayed phase jumps. In fact, phase fluctuations in 2D simulations have been reported to display smaller amplitude than in 1D simulations, where it was suggested that the global movement of the 2D density fluctuations would lead to a smaller apparent change of refractive index than in the 1D case [START_REF] Hacquin | [END_REF].

Despite both the relatively small range of turbulence levels that have been analysed here and the observed scatter in the data, an approximately linear dependence of the phase fluctuations with the turbulence level is recovered for both drift and poloidal plane cut simulations. This indicates that phase jump removal may allow to extend the linear sensitivity of conventional reflectometry phase response, and consequently enable measurements of turbulence levels above what is expected from theory.

Summary and discussion

In this work, the outputs of the gyro-fluid turbulence code GEMR were successfully integrated as inputs to the full-wave code REFMUL. Fixed frequency O-mode simulations were carried out on the LFS mid-plane region allowing to record time traces of synthetic reflectometry I/Q signals. Both the complex amplitude signals and the reflectometry phase signals were retrieved and analysed, for probing frequencies covering the edge and SOL plasma regions where turbulence of ∼3%-6% moderate amplitude develops.

REFMUL simulations with n e data rendered on both drift plane cuts (GEMR coordinate system) and poloidal plane cuts have been performed in order to assess the validity of the usual approximation where the two coordinate systems are not distinguished for practical purposes. Similar trends were observed in both the spectra of the complex amplitude and in the phase fluctuations on the two sets of simulations. For similar turbulence levels, phase fluctuations were higher (>20%) when poloidal planes were used instead of drift planes. Plasma curvature effects may explain this behaviour. For example, in the case of probing wave with frequency f o =30 GHz, an effective radius of curvature r ~0.31 m given by equation predicts the case of strong curvature regime. Considering the sensitivity criterion < k w 2 in the slab geometry and the scaling of equation (4), we can expect reflectometry measurements to have their sensitivity increased from ¢ k 0.4 to ~k 4.4 cm 1 due to curvature effects. One should note that the curvature of the cut-off layer was neglected in the drift planes while the curvature of the incident wave-front was included in both cases. The wave-front curvature alone could according to the same expressions used above. Thus, a significant increase in the sensitivity from that of slab geometry was already expected when considering the 2D drift planes.

The sensitivity could still be further extended by plasma curvature effects, as observed for the poloidal planes. Phase fluctuations arising from the Bragg resonances have been predicted according to a number of models yielding linear or nonlinear relations between the phase fluctuations and the turbulence level at the cut-off location, depending on the amplitude of the density perturbations, e.g. see [23]. Time traces of the reflectometry fluctuating phase components displayed phase jumps of p 2 corresponding to realistic effects observed experimentally. The existence of phase jumps in the simulated data indicate the nonlinear operation regime of the reflectometer. The GEMR turbulence and the REFMUL probing information can be employed to assess the criterion established in equation (5) for the transition between linear and nonlinear regimes. Using the drift plane data, for instance, probing with f o =20 GHz sets the reflectometer operating far from the transition (x = 0.0002) while for f o =32 GHz it is only marginally below (x = 0.98). In fact, probing of GEMR turbulence in drift planes only displayed phase jumps at this probing frequency (as observed in figure 15). However, for the poloidal plane cuts, at least one phase jump occurred in most probing frequencies. Nevertheless, given the small fraction of data representing phase jumps in the whole time series and the marginally nonlinear regime of the reflectometer operation (when the jumps were indeed present), it was chosen to correct them prior to further analysis.

Trends of spectral broadening with increasing probing frequency, in both drift plane and poloidal plane rendered data, were observed. It was possible to retrieve a linear scaling between phase fluctuations and the turbulence levels, as expected in the linear regime of operation, despite the relatively small span of turbulence levels considered here.

For future work, it is envisaged to apply a scaling factor to the density turbulence terms obtained from GEMR in order to probe a larger range of turbulence levels and include values that are closer to what is found experimentally across the SOL, where / dn n e e up to ∼40% is observed [4,13]. This will also allow assessing to which extent it will still be possible to extract undistorted frequency spectra and observe linear phase responses, with and without phase jump removal procedures.

Taking into account the advantages in user/code access and code coupling, for the general use of REFMUL or other reflectometry codes together with GEMR or other turbulence codes, employing a data structure like ITER Integrated Modelling & Analysis Suite (IMAS) should also be considered [77].

Figure 1 .

 1 Figure 1. Cartoon with typical radial profiles of average electron density (solid line) and density fluctuation level (dashed line) at the plasma edge and scrape-off layer (SOL) regions of a medium-sized tokamak. The radial position r LCFS of the Last Closed Flux Surface (LCFS) is used as reference.

Figure 2 .

 2 Figure2. GEMR electron density data on a poloidal plane cut, defined on a polar geometry ( ) q R, , and zoom on region with interpolated data for REFMUL Cartesian grid (x, y).

(

e0(

  background or base plasma), the second set included the density fluctuations ( ) = n x y ,

Figure 3 .

 3 Figure 3. Radial profile from one poloidal plane cut prepared as input for REFMUL. The GEMR full turbulent plasma (solid black) and the base plasma (dashed) are shown with a set of critical density layers corresponding to selected probing frequencies. Also shown are the far-SOL (blue) and adaption layer (red) boundary extensions together with part of the vacuum distance.

Figure 4 .[ 37 ].

 437 Figure 4. Cross-correlation functions (CCF) computed at drift plane cuts. Reference time series were registered at the cut-off locations of probing frequencies f o . An example is shown ( f o =28 GHz) on how to calculate l c directly from the 1/e level.

Figure 5 .

 5 Figure 5. Region of interest of a poloidal plane data matrix from GEMR, together with rectangular spot-size regions centred at the corresponding cut-off locations of probing frequencies f o =28 GHz and f o =32 GHz, of both unfocused (solid) and focused (dashed) antennae.

Figure 6 .

 6 Figure 6. Time traces of average density from drift plane data at selected spot sizes (left) and corresponding power spectra (right).

Figure 7 .

 7 Figure 7. Time traces of average density from drift plane data spotsizes with = w 0.94 cm R

Figure 8 .

 8 Figure 8. Time traces of average density from drift plane data spotsizes with = ẃ m 0.94 cm R R

Figure 9 .

 9 Figure 9. Drift plane radial profiles of average base density (solid black) and of standard deviation of the turbulence level (dotted red). Values at cut-off locations for a set of probing frequencies are highlighted (squares) on each profile.

Figure 10 .

 10 Figure 10. Poloidal plane radial profiles of average base density (solid black) and of standard deviation of the turbulence level (dotted red). Values at cut-off locations for a set of probing frequencies are highlighted (squares) on each profile.

Figure 11 .

 11 Figure 11. Synthetic I/Q signals obtained with drift plane data using f o =28 GHz and f o =32 GHz.

9 .

 9 Results from REFMUL-GEMR simulations Following what was previously prescribed in section 7 regarding the REFMUL simulations, baseband I/Q signals were obtained from sets of REFMUL simulations, at each fixed frequency, using drift planes and poloidal planes separately. Illustrative examples of synthetic [ ] I n and [ ] Q n discretized signals are shown in figure 11.

Figure 12 .

 12 Figure 12. Spectral analysis of the complex amplitude signals ( ) ( ) jA t e i t for a set of probing frequencies f o in drift plane cases.

Figure 13 .

 13 Figure 13. Spectral analysis of the complex amplitude signals ( ) ( ) jA t e i t for f o =29 GHz in poloidal plane cases, with and without a focused antenna.

Figure 14 .

 14 Figure 14. Square root of total integrated power of the complex amplitude signals ( ) ( ) j

Figure 15 .

 15 Figure 15. Synthetic j p and dj phase signals obtained with the corresponding sets of base and fluctuating density matrices of drift plane cuts. 2p plateaus are depicted on the right (broken grey lines).

Figure 16 .

 16 Figure 16. Synthetic phase signals obtained with REFMUL simulations with f o =28 GHz and poloidal plane cuts (left). Correction of dj signal by phase jump removal (right).

Figure 17 .

 17 Figure 17. Fluctuating phase signals using drift plane cuts. Left: synthetic dj signals after phase jump corrections, with imposed offsets of 2π multiples for better visualization. Right: power spectra of the corresponding signals.

( 3 )

 3 can be estimated for the poloidal plane simulations, also yielding that ( )

Figure 18 .

 18 Figure 18. Fluctuating phase signals using poloidal plane cuts. Left: synthetic dj signals after phase jump corrections, with imposed offsets of 2π multiples for better visualization. Right: power spectra of the corresponding signals.

Figure 19 .

 19 Figure 19. Standard deviation of reflectometry phase signals from both planes, against the local turbulence level at the cut-off location. Theoretical values estimated from equation (2) are also shown in dashed lines.

Acknowledgments

This work has been carried out within the framework of the EUROfusion Consortium and of the French Federation for Magnetic Fusion Studies (FR-FCM) and has received funding from the Euratom research and training programme 2014-2018 and 2019-2020 under grant agreement No 633053. The views and opinions expressed herein do not necessarily reflect those of the European Commission. IST activities also received financial support from 'Fundação para a Ciência e Tecnologia' through project UID/FIS/ 50010/2019.