Density wavenumber spectrum measurements, synthetic diagnostic development, and tests of quasilinear turbulence modeling in the core of electron-heated DIII-D H-mode plasmas

Measurements of the turbulent density wavenumber spectrum, δnˆe(k⊥) , using the Doppler Back-Scattering (DBS) diagnostic are reported from DIII-D H-mode plasmas with electron cyclotron heating as the only auxiliary heating method. These electron-heated plasmas have low collisionality, νe∗<1 , Te/Ti>1 , and zero injected torque—a regime expected to be relevant for future fusion devices. We probe density fluctuations in the core (ρ ≈ 0.7) over a broad wavenumber range, 0.5⩽k⊥⩽16 cm−1 ( 0.1⩽k⊥ρs⩽5 ), to characterize plasma instabilities and compare with theoretical predictions. We present a novel synthetic DBS diagnostic to relate the back-scattered power spectrum, Ps(k⊥) —which is directly measured by DBS—to the underlying electron density fluctuation spectrum, δnˆe(k⊥) . The synthetic DBS Ps(k⊥) spectrum is calculated by combining the SCOTTY beam-tracing code with a model δnˆe(k⊥) predicted either analytically or numerically. In this work we use the quasi-linear code Trapped Gyro-Landau Fluid (TGLF) to approximate the δnˆe(k⊥) spectrum. We find that TGLF, using the experimental profiles, is capable of closely reproducing the DBS measurements. Both the DBS measurements and the TGLF-DBS synthetic diagnostic show a wavenumber spectrum with variable decay. The measurements show weak decay (k −0.6) for k < 3.5 cm−1, with k −2.6 at intermediate-k ( 3.5⩽k⩽8.5 cm−1), and rapid decay (k −9.4) for k > 8.5 cm−1. Scans of physics parameters using TGLF suggest that the normalized ∇Te scale-length, R/LTe , is an important factor for distinguishing microturbulence regimes in these plasmas. A combination of DBS observations and TGLF simulations indicate that fluctuations remain peaked at ITG-scales (low k) while R/LTe -driven TEM/ETG-type modes (intermediate/high k) are marginally sub-dominant.


Introduction
In tokamak confinement, understanding plasma fluctuations responsible for the transport of heat and particles across flux surfaces remains a fundamental hurdle for fusion energy.Steady-state kinetic profiles, i.e. density and temperature (n, T respectively), determine the performance of plasma discharges.The profiles themselves are the result of balancing sources (e.g.auxiliary heating, gas puffing) and losses (e.g.radiation, edge recombination) mediated by transport processes.Decades of research has shown that significant transport is caused by plasma micro-instabilities [1].Identifying unstable modes and validating reduced models with fluctuation diagnostics is a vital step for building our understanding of plasma performance and making predictions.
In this work we use the Doppler Back-Scattering (DBS) diagnostic to investigate the character of density fluctuations in the core of H-mode plasmas with dominant electron heating and zero injected torque.During the steady-state phase of these H-mode discharges, the only auxiliary heating 5 is from electron cyclotron resonance heating (ECH).This plasma regime is expected to be relevant to future burning fusion devices where fusion-born alphas primarily heat electrons and Neutral Beam Injection (NBI) applies a smaller torque (relative to today's tokamaks).We present a unique experimental dataset from scans of the DBS beam's probing angle across repeated discharges.The primary experimental result is the back-scattered power, P s , over a broad range of wavenumbers perpendicular to the background magnetic field, k ⊥ .
The level of back-scattered power at each measured wavenumber is a function of the fluctuating density wavenumber amplitude spectrum, δn(k ⊥ ).The δn(k ⊥ ) spectrum is an important quantity for diagnosing microinstabilities and indirectly assessing their role in transport.Furthermore, δn(k ⊥ ) can be predicted by theory or simulations; providing an avenue to compare with measurements.
We use the beam-DBS model developed by Hall-Chen et al [2] to explore the relationship between the directlymeasured P s (k ⊥ ) spectrum and the underlying δn(k ⊥ ).We use the SCOTTY beam-tracing code to calculate the DBS weighting function along the probing beam trajectory.The weighting function is used in our novel synthetic DBS diagnostic which takes an analytic or code-based approximation of δn(k ⊥ ) as input and predicts a synthetic back-scattered power, P syn.s .We also invert this calculation, such that the measured P s is an input, and the parameters of a model δn(k ⊥ ) are fit iteratively.The final product of this inverse calculation is an estimate of the underlying δn(k ⊥ ) spectrum which reproduces the measured P s (k ⊥ ).
For the plasma presented in this work, we use the Trapped Gyro-Landau Fluid (TGLF) code to approximate the saturated δn(k ⊥ ) spectrum.Because TGLF is a reduced-model quasilinear fluid code, we take an additional step to compare the eigenvalue spectrum from TGLF with the linear eigenvalue spectrum from the gyrokinetic code CGYRO.
Leveraging the computational speed of beam-tracing and the TGLF reduced model, we perform scans of physics parameters.These scans help us identify parameters which contribute significantly to the synthetic diagnostic output, i.e. parameters which are predicted to affect the observed P s spectrum.We find that small changes (±10%) in the inverse-∇T e scale length, R/L Te , cause a significant change in the TGLF-δn(k ⊥ ) spectrum.This change in the predicted fluctuation spectrum carries through to the P syn.s spectrum.The combination of our measurements and synthetic diagnostic modeling suggests R/L Te -driven modes (TEM/ETG) are marginally subdominant to ion-scale modes (ITG) in this plasma.
This paper is organized as follows: section 2 describes the theory of measuring δn with DBS including the development of the synthetic DBS diagnostic.Section 3 describes the codes used for modeling and the method for estimating δn using TGLF.Section 4 discusses the experiment and section 5 focuses on the specifics of data analysis.Sections 6 and 7 contain the results and analysis/discussion respectively.Finally, appendix A provides an abridged dictionary of symbols used in this paper.

DBS background and theory
DBS is a fusion plasma diagnostic capable of measuring turbulent flow velocities [3,4], density fluctuation spectra [5], and radial correlation lengths [6].DBS is an established diagnostic used on tokamaks and stellarators around the world [4,[7][8][9][10].This section will focus on the theory of measuring the density fluctuation wavenumber spectrum by varying the launch-angle of the DBS probing beam.
To understand the theoretical DBS response to varying the launch-angle we begin with a description of the DBS scattering geometry.Then, we use a simplified model of DBS adapted from Ruiz et al [11] to illustrate the underlying physics allowing for a measurement of the fluctuation wavenumber spectrum.Later we will introduce and use the beam model of DBS developed by Hall-Chen et al [2].Both models assume that density fluctuations are sufficiently weak that diagnostic response is linear, see appendix B for more regarding the linearity assumption.

Scattering geometry
We focus on perpendicular wavenumbers, k ⊥ , motivated by the fact that instabilities in tokamaks tend to have strong anisotropy with the respect to the background magnetic field such that turbulent structures satisfy k ∥ ≪ k ⊥ .In the configuration used here for DBS, the scattering is dominantly due to the short perpendicular scale of the turbulence and the parallel structure does not contribute significantly.Therefore, in our modeling we ignore parallel structure and focus on the 2D plane perpendicular to the magnetic field.
Figure 1 illustrates various coordinate systems common in tokamaks along with a representation of DBS scattering geometry.Figure 1(a) illustrates the intersection of a probing DBS ray (red) with a flux surface in 3D.In figure 1(a) we see the magnetic field, B, piercing a plane labeled (b).A projection of this plane is expanded as figure 1(b).Note that the flux surface in figure 1(a) and the expanded regions in figures 1(b) and (c) are in the vicinity of the turning point of the DBS trajectory.
The plane perpendicular to B in figure 1(b) defines the normal ên , and binormal êb directions.Anywhere in the plasma, wavevectors perpendicular to the magnetic field can be written, 1(b) and (c) illustrate the scattering processing using perpendicular (purely binormal) wavevectors at the turning point.The wavevector of the incident DBS ray, k i scatters from density fluctuations with k ñ to produce an outgoing ray with wavevector k s .In this situation, the measured wavevector, k meas.= k ñ = −2k i by the Bragg scattering relation.This vector relation is a general expression and is not limited to the vicinity of the cutoff.
Figure 1(c) is provided to illustrate the scattering vectors when projected onto an unwrapped flux-surface.It is important not to confuse the binormal direction (subscript b) with the direction along the magnetic field-line (parallel symbol, ∥).The pitch-angle of the magnetic field line is exaggerated to differentiate the binormal (ê b ) and poloidal (θ) directions.Only in the limit of vanishing magnetic pitch would the binormal direction coincide with the poloidal direction.

Simplified model DBS
The back-scattered wave amplitude, A s can be expressed as a real-space integral of the DBS weighting function U, and the fluctuating density δn [11,12], Within this simplified model we treat the measured wavevector, k meas., as a constant.The domain of the integral in equation ( 1) is essentially set by U(r) which follows the trajectory of the probing DBS beam.The scattered power is calculated via We then perform an ensemble average, exploiting many realizations of the fluctuating density.We assume δn is a stationary random process such that the autocorrelation function, R, depends only on the difference in position, i.e. ⟨δn(r)δn * (r ′ )⟩ = R(r ′ − r). 6The autocorrelation function is related to the spectral density function, S, by Fourier transform, To write the average back-scattered power in k-space we substitute equation (3) for ⟨δn(r)δn * (r ′ )⟩ in equation ( 2) to find, Exchanging the order of integration we can group real-space terms (involving r, r ′ ) to write, We recognize the term inside the modulus-square to be the Fourier transform of the weighting function F[U] = Û (we use the 'hat' symbol,ˆ, to denote a Fourier amplitude).Importantly, the Fourier transform of the weighting function is evaluated at k − k meas., Finally, one can show that the stationary random process assumption implies that the spectral density function is equivalent to the averaged, modulus-square of the (Fourier-space) fluctuating density, S(k) = ⟨|δn(k)| 2 ⟩.The δn(k) spectrum is written as a function of the vector k to allow for anisotropy in the perpendicular (k n , k b ) plane.In earlier expressions we used δn(k ⊥ ) for brevity; the more precise expression is δn(k) = δn(k n , k b ) (recall that we ignore the parallel direction in this subsection).Equation (6) states that the back-scattered power as a function of the measured wavenumber, ⟨P s (k meas.)⟩, is essentially a convolution of the fluctuating wavenumber spectrum, ⟨|δn(k)| 2 ⟩ with the DBS k-space weighting function, Û.When we scan the launch-angle of the DBS beam, we change the measured wavenumber k meas. .Thus, a launch-angle scan systematically translates the DBS weighting function in k-space to probe different regions of ⟨|δn(k)| 2 ⟩.In the (non-physical) limit where Û(k − k meas. ) ∝ δ(k − k meas. ) we would obtain a 'perfect' measurement of ⟨|δn(k meas.)| 2 ⟩.
The assumptions built into the simplified model from the previous section miss important physical effects which can impact the measurement of the density wavenumber spectrum.To more accurately model DBS measurements, we leverage the beam-DBS theory developed by Hall-Chen et al [2].
In beam-DBS theory we treat the probing field as a Gaussian beam and we use the reciprocity framework [12] to derive an expression similar to equation ( 6), but with a more complicated DBS weighting function U related to how well the probing beam illuminates a scattering volume in the plasma.However, the essential physics of a k-space convolution between Û and ⟨|δn(k)| 2 ⟩ remains relevant.
With a Guassian probing beam, the measured wavevector k meas. is not fixed; rather it evolves along the trajectory of the beam.Furthermore, the beam-model includes diagnostic effects such as wavevector 'mismatch' where poor alignment with the binormal direction attenuates the scattered signal.This effect was studied on DIII-D and validated against the beam model in [13,14].
The average back-scattered power, using the beam-DBS model, can be written as an integral along the trajectory of the probing beam, parameterized by the path-length l, cf equation (196) in [2].In accordance with the original work we label individual terms in the integrand of equation ( 7): (1) U p : the 'polarization' term, (2) U b : the 'beam' term, (3) ∝ g −2 : the 'ray' term, (4) Gaussian in θ m : the 'mismatch' term, and (5) ⟨|δn| 2 ⟩: the 'turbulence' term.The derivation and complete definition of each term can be found in the original text [2].The integral in equation ( 7) is taken along the trajectory of the beam through the plasma with all terms evaluated at the measured wavevector k meas.(l) satisfying the Bragg scattering condition at each point along the beam path.The measured wavevector points in a direction (labeled ê1 ) perpendicular to the background magnetic field and coplanar with the beam group velocity, g.To properly evaluate the density fluctuation spectrum along the beam path, we perform an additional vector decomposition of k meas.(l) into its normal and binormal components, Wherein k n,meas.and k b,meas.are the measured normal and binormal wavenumbers, respectively.Note that every term in equation (8), including the unit vectors, varies along the path of the beam trajectory and can be parameterized by l.Generally speaking, DBS measurements are localized in the vicinity of cutoff surfaces.Although commonly assumed a priori for DBS analysis, this statement will be justified later with results from the plasma studied here.At cutoff (l = l c ), the normal component of the beam vanishes such that k meas.ê1 | lc = k b,meas.êb | lc .Throughout this work, when a single value of k meas. is assigned to a given trajectory, it is implied that k meas. is the measured, binormal wavenumber at cutoff.Importantly, it is not appropriate to assign a simple, local uncertainty on the measured wavenumber (∆k meas. ) when using the integral in equation (7).Localized ∆k meas.values cannot be separated from real-space integration along the beam trajectory (for more discussion see section 9.1 of [2]).
The total DBS weighting function is a combination of four terms, U ∝ U p U b g −2 exp (−2θ 2 m /∆θ 2 m ).The polarization term, U p , takes different forms depending on whether the probing DBS beam is polarized O/X-mode.For X-mode, U p depends on the relative sizes of the electron cyclotron frequency and the electron plasma frequency whereas for O-mode, U p is roughly constant.
The mismatch term, exp (−2θ 2 m /∆θ 2 m ), also varies along the beam trajectory, but is primarily set by the initial launch angle(s) relative to the plasma magnetic equilibrium.The mismatch angle is given by θ m ≡ sin −1 ( k • ê∥ ).Improper aiming of the probing DBS beam results in a Gaussian attenuation of the back-scattered signal.When mismatch is not accounted for, measurements of the scattered power can be compromised (i.e.much lower than the real value).The parameter ∆θ m sets the width of the Gaussian mismatch attenuation.Generally speaking, ∆θ m is smaller for larger k meas., i.e. when making higher-k meas.measurements it becomes increasingly important to have θ m < ∆θ m at cutoff.For the cases presented here, In cases where θ m ≈ 0, The ray term, ∝ g −2 , varies along the beam trajectory as k −2 .Thus, when the beam reaches cutoff and the normal component of k vanishes, k −2 is maximized.This effect creates an Airy-type enhancement of the probing field at cutoff, contributing to the localization of the measurement.This k −2 component of the weighting function is a central part of established DBS measurement theory [15].
We neglect the 'beam' term, U b , in the integrand of equation (7) with the assumption that the beam properties do not vary significantly while scanning the DBS launch angles.Properly capturing the behavior of the beam near cutoff for a wide array of launch-angles is an ongoing area of research [16].
To compute terms in the integrand of equation ( 7) (except the turbulence term), and to evaluate the 3D beam trajectory, we use the beam-tracing code SCOTTY [2].This code was developed as part of the beam-DBS theory and internally calculates components of the DBS weighting function.

Synthetic DBS diagnostic
Systematically scanning the DBS launch-angle directly measures the back-scattered power wavenumber spectrum, P s (k meas.).However, this is an indirect measurement of the underlying fluctuating density wavenumber spectrum, δn(k).To bridge the gap between P s (k meas. ) and δn(k), we rely on synthetic diagnostic modeling.We apply our synthetic diagnostic to solve both the forward and inverse problems: On the one hand, the forward problem is a straightforward evaluation of equation ( 7) provided one has access to δn(k).However, modeling δn(k) via simulations can be prohibitively expensive for two reasons: (1) the saturated fluctuation spectrum is an inherently nonlinear quantity, and (2)-depending on the range of probed wavenumbers-multi-scale runs may be necessary for a complete comparison with measurements.
To forgo multi-scale nonlinear simulations, we present a novel method to obtain an approximate δn(k) using the reduced model TGLF.This will be covered in section 3.
On the other hand, the inverse problem is solved by first assuming a particular, analytic model for δn(k; ⃗ µ) where ⃗ µ represents a set of model parameters.Then, we iterate ⃗ µ until the model δn reproduces the experimental P s (k meas. ) via equation (7).This approach introduces bias through the functional form of the model δn and may converge for several different models.Convergence indicates that the δn spectrum, with optimized parameters ⃗ µ * , can reproduce the observations.

Quasi-linear modeling
In reality, the fluctuating density wavenumber spectrum, δn(k), is a result of the saturated, nonlinear state of unstable plasma modes.When it comes to approximating δn(k) numerically, one cannot calculate δn(k) with linear theory alone.Unfortunately, performing global multi-scale nonlinear gyrokinetic simulations is prohibitively expensive for routine analysis.Instead, we use the quasi-linear TGLF [17] model to compare our DBS measurements with theoretical predictions.TGLF is a fluid code solving a linear eigenvalue problem over a grid of binormal wavenumbers, k y , to determine the frequency, ω, and growth rate, γ, of unstable plasma modes.To calculate fluxes and fields such as the fluctuating density, potential, and temperature; TGLF uses a saturation rule informed by a database of nonlinear gyrokinetic simulations [18].In essence, the saturation rule models the fluctuating potential spectrum and the quasi-linear weights (phase shifts) necessary for calculating transport fluxes.
In this work, we use the SAT2 rule of TGLF [19].The SAT2 model improves upon previous rules (e.g.SAT0, SAT1) with a 3D model for the saturated potential spectrum, δ φ(k x , k y ; θ), wherein k x is the radial/normal wavenumber, k y is the binormal wavenumber, and θ is the poloidal angle.Because TGLF is a quasi-linear model, the most unstable linear modes at each k y are used to approximate the saturated nonlinear state of the plasma.To verify, for our case, that TGLF is able to accurately calculate the linear eigenvalue spectrum, we will compare with the linear spectrum calculated by the gyrokinetic code CGYRO [20].Generally speaking, the latest available saturation rule is best equipped to mimic the higher fidelity gyrokinetics.
When experimental data is used to test simulations, any differences in definition between the theoretical and measured wavenumbers must be addressed.Both TGLF and CGYRO use the same wavenumber definitions and normalizations [20]: where ρ s,unit ≡ c s /Ω D,unit is the effective ion-sound gyroradius, c s = T e /m D the ion sound speed, Ω D,unit = eB unit /(m D c) is the ion cyclotron frequency based on an effective magnetic field strength B unit ≡ (q/r)dψ/dr (m D is the mass of Deuterium, ψ is the poloidal flux).The k x wavenumber is radial with mode number p and L representing the radial domain size (e.g. a flux tube in CGYRO).The k y is binormal with toroidal mode number n, flux surface minor radius r, and safety factor q.
Ruiz et al [11] provides the necessary mappings between these theoretical (k x , k y ) wavenumbers and normal, binormal wavenumbers (k n , k b ) which correspond to physical inverselengths.The relationship for up-down symmetric (UDS) flux surfaces at the outboard midplane (θ = 0) is given by, where R is the major radius of the flux surface, α is the Clebsch coordinate of the magnetic field 7 , and κ is the flux surface elongation.For a more detailed discussion about the wavenumber mappings see appendix C. When applied, equation ( 10) scales the theoretical radial and binormal wavenumbers (i.e.k x , and k y ) to account for flux-surface shaping effects.In the radial direction, Shafranov shift compresses flux surfaces (k n > k x /ρ s,unit ) whereas plasma elongation expands the binormal direction (k b < k y /ρ s,unit ).
To address the forward problem described in section 2.4 we leverage TGLF's 3D spectral-shift model for the saturated potential spectrum, (11) cf equation ( 34) of [18].This model is used by TGLF to approximate the saturated potential spectrum across a database of nonlinear gyrokinetic simulations.This particular equation is an extension by Staebler et al [18] to include the effects of E × B shear, γ E×B into a previous model for δ φ(k x , k y ; θ) presented in [19].The effects of γ E×B appears, in part, through the k x -shift parameter, k x0 .In equation (11), G(θ) is the only θdependent term.To simplify, we take θ = 0, focusing ourselves to the outboard midplane (where DBS measurements are primarily made, and where equation ( 10) is applicable).
The numerator of equation (11), δ φ(k y )| kx=0,θ=0 , is the peak of the potential spectrum in the absence of equilibrium γ E×B .The value of δ φ(k y ) is calculated using TGLF model parameters and the growth rates, γ ky , of the dominant, linearly unstable modes.The denominator of equation ( 11) creates a Lorentzian in k x .The Lorentzian k x spectral shape was shown to be a good fit to the radial wavenumber spectrum from nonlinear gyrokinetic simulations [21].The properties of the k x -Lorentzian: k model x , k x0 are functions of k y .The terms α x , and σ x are constants of the saturation rule.The precise definitions of all parameters in equation ( 11) can be found in [18].
We use equation (11) to reconstruct the 2D (k x , k y ) saturated potential spectrum, δ φ(k x , k y ; θ = 0) from the outputs of TGLF.The quasi-linear weights, W QL ky are used to connect the potential spectrum and the density spectrum with, This assumes that the QL weights themselves are not functions of k x .The quasi-linear weights are defined by, wherein the numerator and denominator are the 1D gyro-Bohm normalized intensity spectra outputs from TGLF (cf equation ( 4) of [18]).Lastly, we convert from (k x , k y ) to (k n , k b ) via equation ( 10) yielding an approximation of the 2D saturated δn 2 (k n , k b ) spectrum.
For use in the synthetic DBS diagnostic, the approximate 2D saturated δn 2 (k n , k b ) spectrum is inserted into the integrand of equation (7) where it is evaluated using the measured normal, binormal wavenumber components (equation ( 8)).

Experiment
DBS measurements of the fluctuating density wavenumber spectrum were made in DIII-D ECH H-mode discharges (Hmodes with ECH as the only auxiliary heating).This regime of operation is relevant for studying plasmas with dominant electron-heating and zero injected torque.In future, burning reactors the fusion-born alphas will predominantly heat electrons.These future fusion reactors are also expected to have reduced reliance on NBI, and are therefore expected to have low to zero injected torque (the extent to which intrinsic torque, especially in RF-heated plasmas, is able to produce rapid rotation is the subject of ongoing research [22,23]).
Figure 2 illustrates the evolution of a typical ECH H-mode discharge from this study (DIII-D 189998).The plasma current and magnetic field are held constant at I p = 0.9 MA, and B ϕ = −1.8T (i.e.opposite the I p direction).The edge safety factor, q 95 ≈ 6, and the normalized plasma beta, β n ≈ 0.9, produce no significant MHD activity.The normalized confinement factor H 98,y2 ≈ 1.These are diverted discharges in the lower-single-null (LSN) configuration with plasma elongation, κ = 1.7.
Figure 2(a) shows the line-average density, ⟨n⟩, and the plasma current, I p . Figure 2  Once the I p -flattop is established, 2.7 MW of NBI is used over 400 ms to assist the LH transition (see figures 2(a) and (b)).This NBI pulse uses two sources in opposing (tangential) directions such that near zero torque is applied to the plasma.Halfway through the NBI pulse, 2.1 MW of ECH was applied near the core of the plasma (ρ ∼ 0.3, see figure 3(b)) and held constant for the remainder of the discharge.This sequence of external heating was used to trigger a reproducible H-mode transition while imparting near zero torque on the plasma.We were reliably able to sustain the H-mode with ECH as the only external heating source.
Neutral beam 'blips' were used to obtain data from beamdependent diagnostics such as Motional Stark Effect (MSE) and Charge Exchange Recombination Spectroscopy (CER).The blip from a single NBI source lasts 10 ms.Each pulse visible in figure 2(b) consists of two simultaneous blips (both are directed co-I p ). Blips are spaced 400 ms apart, resulting in only sparse data from both MSE and CER.The long interblip time minimizes the influence of heating, torque, and fast particles from NBI.
These discharges do not use gas feedback to control the line average density.The density rises to its maximum value after the LH transition.Then, over the next 400 ms, the density drops due to a combination of global ECH density pump-out, high ELM frequency, and the lack of a particle-source from NBI.The plasma finally enters a quasi-stationary phase with lower ELM frequency and a density of ⟨n⟩ ∼ 2.8 × 10 19 m −3 near 2500 ms.
Figure 3 shows an array of representative kinetic profiles and their inverse gradient scale-lengths, a/L X ≡ −(a/X)∇X, where a is the minor radius of the last closed flux surface.The DIII-D Thomson scattering (TS) system and the Profile Reflectometry system (Refl.)[24] are used in combination to fit the n e profile (figure 3(a)).Similarly, the electron-cyclotron emission diagnostic [25] is used in combination with TS to fit T e (figure 3(b)).Impurity-CER (impCER) provides measurements of the carbon impurities in the plasma (figure 3(d)) while main-ion CER (miCER) measures the main-ion (deuterium) T i profile (figure 3(c)) [26].The effective charge-state of the plasma, Z eff ≈ 2. Profile fits are performed using the OMFITprofiles module [27].The uncertainties in the fits are calculated with linear error propagation theory, accounting for the correlation between fit-model parameters.
Due to high electron-heating and relatively low lineaveraged density, the plasma achieves very low electron collisionality.In the core ν * e < 0.1, while the pedestal has ν * e < 1.We define ν * e ≡ ν i,e (qRv −1 e ϵ −3/2 ), where ν i,e is the electronion collision frequency, q is the safety factor, ve is the electron thermal velocity, ϵ is the aspect ratio, and R is the flux surface major radius.In terms of kinetic quantities, ν * e ∝ n e T −2 e .The plasma also possesses a large ∇T e , yielding a gradient-ratio value of η e ≈ 3 in the radial regime of interest (ρ ∼ 0.7).The gradient-ratio is defined as, η s ≡ L ns /L Ts = (n s ∇T s )/(T s ∇n s ).Lastly, the core T e /T i peaks at 2.0, and remains >1 everywhere inside the pedestal-top.Without NBI, the main ions are heated by collisions with electrons.The electron-ion heat source is largest in the core (where T e and T i differ the most) and extends to the pedestal-top.

DBS launch-angle scans
DIII-D has two eight-channel, V-band DBS systems with fixed frequencies ranging between 55-75 GHz [28].These systems are toroidally separated on the machine-one at the 240 • port and the other at the 60 • port.We will refer to these independent systems as DBS240, and DBS60.The toroidal location of the systems is not relevant other than as a label.Importantly, the DBS240 system possesses the ability to steer the probe-beam both toroidally and poloidally [29].The toroidal steering capability of the DBS240 system allows us to match the magnetic pitch-angle and avoid mismatch attenuation (discussed previously in section 2.3) at higher-k [14].The DBS60 system has a fixed toroidal launch-angle of 0 • making it suitable for lower-k measurements.
The launch-angles of both DBS systems remain constant during a plasma discharge.Therefore, we repeat discharges to scan the launch-angles (and therefore the measured wavenumber).For the wavenumber scan presented here, we repeated the same ECH H-mode discharge a total of six times: 189993,189994,189997,189998,190129,190315.All repeats are in good agreement with each other with the exception of 190129, which back-transitioned to L-mode before the I p ramp-down.To build a cohesive wavenumber spectrum by combining measurements from both DBS systems, one discharge (189998) was used to cross-calibrate the systems by targeting the same wavenumber.To validate this technique, we matched the probed wavenumber of both systems a second time (190129) to ensure that the calibration factor derived at one wavenumber was still applicable at another.

Data analysis
The analysis of DBS data consists of (1) frequency-integration of the Doppler-shifted component of the scattered signal power spectral density and (2) ray/beam tracing to identify the measured wavenumber and cutoff location.Together these provide the back-scattered power wavenumber spectrum, P s (k meas.).
All DBS data presented here is from the same frequency (72.5 GHz) channel of each multichannel system.Using a single channel allows us to avoid cross-channel calibrations and more directly compare data between discharges.As mentioned in section 4, specific discharges were used to address cross-system calibration, necessary to compare the 72.5 GHz channels of the independent DBS60/240 systems.

Spectral analysis
For each discharge in the wavenumber scan, DBS data is sliced into multiple (generally 3-4) 150 ms windows centered 300 ms away from NBI blips.These time windows were selected to minimize the influence of NBI.The digitized (5 MHz sampling) quadrature mixer outputs, I(t) and Q(t), are cast into a complex-valued representation of the scattering signal, S(t) = I(t) + iQ(t).The signal S(t) is processed with a shorttime Fourier transform using 16 384 points, Hanning windows, and 50% overlap between realizations to produce a spectrogram.The spectrogram is then smoothed in the frequency domain by a moving-average with a 10 kHz window.
Figure 4(a) provides a representative 72.5 GHz DBS spectrogram from discharge 190129.Figure 4(c) overplots three time-averaged spectrograms from discharges 190129, 189997, and 189994.We observe the Doppler shift increasing and the spectral power decreasing as we increase the measured wavenumber.This increase in the observed Doppler shift is expected as the Doppler shift is proportional to the measured wavenumber, i.e. in general ω dopp.= k meas.• v.
Beyond minimizing the influence of NBI blips with our choice of 150 ms time windows, we also remove the influence of ELMs on the DBS scattered-power calculation by performing an ELM phase-locked analysis.Every 150 ms window of DBS data is synchronized to the corresponding D α filterscope signal (see figure 4(b)).We retain time-slices of each spectrogram corresponding to 50%-95% ELM-phase only.These time-slices are shown as shaded regions in figure 4(b).The ELM-filtering process reduces the total number of spectrogram time-slices by approximately 50%.The remaining 'valid', i.e. not influenced by ELMs, time-slices are analyzed to identify the Doppler-shifted component and its spectral power.Where necessary, the Doppler-shifted component is isolated by applying a high-pass notch filter to remove spectral power near zero frequency.Finally, spectra are fit with a Gaussian and the spectral power is calculated by numerical integration.

Ray tracing
We performed 3D ray tracing simulations to determine the measured wavenumber, k meas., for each 150 ms window of DBS data (generally 3-4 windows per discharge).Ray tracing simulations are performed with the GENRAY code [30] provided with magnetic equilibrium reconstruction from the EFIT code [31] and fits to electron density profile data.
X-mode polarization is used exclusively for the wavenumber scan discharges.Using consistent polarization avoids additional corrections and provides for a simpler comparison across discharges.The subtleties of combining O/X-mode wavenumber scan data was explored extensively by Happel et al [5] with related modeling by Lechte et al [32].
Figure 5 illustrates the ray tracing results for the wavenumber scan.Figure 5(a) provides an annotated poloidal cross section of the plasma with a subset of the ray trajectories from the wavenumber scan.Notably, we show the lowest probed wavenumber and the highest probed wavenumber along with the matched-wavenumber case used for crosssystem calibration.
Figures 5(a) and (c) show the tendency of higher launch angles to probe larger radii.This is an expected effect due to refraction and increasingly oblique incidence relative to the X-mode cutoff.This effect broadens distribution of probed-ρ locations despite constant launch frequency (72.5 GHz).Furthermore, small variations in the density profile across discharges and over time within any one discharge, also broadens probed-ρ.Figure 5(d) histogram of the probed locations for the entire wavenumber-scan dataset.The mean and standard deviation of the probed radial locations are ρ 0.67 ± 0.07.The mean probed location is used later for more extensive analysis.

Results
The primary experimental result is the back-scattered power wavenumber spectrum, P s (k meas. ) shown in figure 6.We observe a spectrum with variable spectral decay: for lower wavenumbers (k < 3.5 cm −1 ) the spectrum is nearly flat with P s ∼ k −0.6 , for intermediate wavenumbers (3.5 < k < 8.5 cm −1 ) the spectrum decays with P s ∼ k −2.6 , and for high wavenumbers (k > 8.5 cm −1 ) we observe steep decay with P s ∼ k −9. 4 .
The variable slope in our P s (k meas. ) spectrum motivated a nonlinear fit using a generalized exponential function, The optimal P s (k meas. ) fit has parameters: b = 0.8, c = 1.9, and d = 1.1.Interestingly, the proximity of d to 1.0 implies an almost purely exponential scattered-power spectrum.The secondary x-axis above figure 6 (also used in later figures) shows the mapping between k meas.and the theoretical k y value.Equation (10a) is used to calculate k y given the measured binormal wavenumber, k b = k meas. .The plasma parameters in equation (10a) (κ, q, ρ s,unit etc) are evaluated using local quantities for each data point in figure 6.Thus, this mapping includes the radial variation of said plasma parameters over the entire dataset.The overall k y vs. k meas.trend for this plasma was found to follow k y ≈ 0.44(k meas. ) 0.92 .
Similarly, to normalize the k meas.wavenumber by the experimental ρ s , we use the local ρ s value corresponding to each cutoff location.This provided an overall k meas.ρ s vs. k meas.trend of k meas.ρ s ≈ 0.45(k meas. ) 0.84 .Purely by coincidence, when the trends for k y and k meas.ρ s are combined, we find k y ≈ 1.1(k meas.ρ s ) 1.1 .Therefore, the upper k y axes on figure 6, and subsequent figures, can roughly be used as k meas.ρ s .It is important to note that generally k meas.ρ s ̸ = k y .
The DBS data in figure 6 is re-plotted against k meas.ρ s in figure 7 using a linear x-axis.In figure 7 we also observe variable decay.For k meas.ρ s < 3.5 the spectrum decays exponentially with ζ ≈ 1.7 (ζ : P s ∼ e −ζ(kmeas.ρs)).For k meas.ρ s > 3.5 the spectrum decays with ζ ≈ 3.9.Once again, the theory k y value is provided as a secondary x-axis above figure 7.
The experimental P s (k meas. ) spectra shown in figures 6 and 7 are indirect measurements of the underlying fluctuating wavenumber spectrum, δn(k n , k b ).Additional modeling to address the forward and inverse problems introduced in section 2.4 is presented next.All detailed analysis in the following subsections uses profiles and the magnetic equilibrium from a single, representative discharge/time: DIII-D 189998, 3000 ms.Back-scattered power wavenumber spectrum.The DBS data follows the same marker/color scheme as figure 5(c).Black dashed lines represent linear (in log-log) fits with each spectral index provided by annotations.The linear fits, which normally pass through the data points, have been displaced to make them visible.The red line is the result of a nonlinear fit using the function in equation ( 14).Note the secondary x-axis above the plot providing a mapping between the real-space kmeas., and the normalized ky.

Forward-problem results
In this subsection we show the results of forward-modeling the DBS back-scattered power wavenumber spectrum.Recall, one has the freedom to provide any model δn(k n , k b ) spectrum and compute a synthetic back-scattered power wavenumber spectrum, P syn.s (k meas.).In this work we use the TGLF code to to approximate the δn(k n , k b ) spectrum.We use the symbol P TGLF s (k meas. ) to denote the synthetic DBS back-scattered power wavenumber spectrum based on TGLF.Given that DBS measurements come from a spread of radial locations (see figure 5(d)), we include the radial variation in the TGLF-δn spectrum in our forward-model.Motivated by the assumption that DBS signals are localized to cutoff surfaces (this assumption is verified later), we first identify the cutoff location for a given DBS trajectory (ρ = ρ c ), then we run TGLF at ρ c , and use the TGLF δn(k n , k b ; ρ c ) to calculate P TGLF s .For a detailed look at the TGLF results we will focus on the averaged probed location (ρ = 0.67) and return to the radial variation later.To verify the eigenvalue spectrum from TGLF is in agreement with linear gyrokinetics, we perform a scan of k y with CGYRO (at ρ = 0.67).Figure 8 shows the growth rates, γ and real frequencies, ω from TGLF and CGYRO.
Generally we find good agreement between the two codes.Both codes agree that for k y < 0.7, modes propagate in the iondrift direction with growth rates, γ ≲ 0.1(c s /a).For k y > 0.7 modes propagate in the electron-drift direction with increasing γ.The intermediate-scale region, just above the transition from the ion-direction to the electron-direction, is also where the TGLF and CGYRO disagree the most in terms of γ.
Figure 9 provides an image of the TGLF-δn 2 (k n , k b ) with log-scaled color at ρ = 0.67.As the DBS beam propagates from the edge plasma to the cutoff, then back outwards, the beam samples different turbulent wavenumbers.The line superimposed on the image in figure 9 is a projection of a single DBS k-space trajectory.If we superimpose a different DBS trajectory with different launch angles, the line in figure 9 would translate horizontally.Thus, varying the DBS launch angles allows us to sample roughly vertical (approximately fixed k b ) slices of the δn spectrum.The location where the DBS trajectory intersects the k n = 0 axis sets the measured binormal wavenumber at cutoff, k meas. .
Figure 10 shows several components of the DBS weighting function calculated using the SCOTTY code (see equation (7)).The case shown in figure 10 corresponds to the k-space trajectory shown in figure 9.In particular, we show the 'ray' (figure 10(a)), 'polarization' (figure 10(b)), and 'turbulence' (figure 10(c)) terms of integrand in equation (7).Recall that the ray component is ∼ k −2 (when θ m ≈ 0; with k being the wavenumber of the central ray) and provides some localization for the signal.
The localization provided by the spectrum-shown in figure 10(c)-can be explained by the anisotropy of the δn 2 spectrum in figure 9.The δn 2 spectrum generally exhibits low amplitude at higher radial wavenumbers.The absence of large-|k n | power in the δn 2 (k n , k b ) spectrum plays a significant role in localizing DBS measurements.
The (k n , k b ) anisotropy of the δn spectrum along with the DBS weighting function support the approximation that the majority of the scattered power originates from a relatively narrow region near the cutoff.Thus, radial variation in the TGLF-δn spectrum is included by running TGLF at various radial locations corresponding to the cutoff of each DBS trajectory.
The complete P TGLF s (k) spectrum is calculated by first running a series of beam-tracing simulations (representing the launch-angle scan).Separately, we reconstruct the TGLFδn(k n , k b ) spectrum from the TGLF simulation performed at the cutoff location.Lastly, we calculate the integral in equation ( 7) over the path of the beam.Each beam tracing simulation provides one k meas.and allows us to calculate one P TGLF s value.When these steps are repeated over the entire set of beam-tracing simulations we produce the P TGLF s (k meas. ) curve.The final P TGLF s (k meas. ) spectrum is shown in figure 11.
To treat the density fluctuations consistently over the launch-angle scan, we do not normalize each P TGLF s separately.However, an overall normalization is applied such that the peak P TGLF s = 1.Due to this final multiplicative P snormalization, the amplitude of the P TGLF s (k meas. ) curve has arbitrary units (a.u.).Given the back-scattered power in absolute units at one wavenumber, the curve would be translated vertically (in log-space).Importantly, the relative magnitude of neighboring wavenumber points is fixed and the radial variation in δn(k n , k b ), as predicted by TGLF, is included.This allows us to quantitatively compare the shape of the predicted P TGLF s (k meas. ) spectrum with measurements.

Inverse-problem results
To solve the inverse problem outlined in section 2.4, we use the experimental P s (k meas. ) spectrum as the target of a nonlinear least-squares optimization.We use the symbol P inv.
s (k meas. ) to represent the synthetic back-scattered power s (k) used equation ( 15) for the underlying δn spectrum.The DBS data shown in figure 11 is identical to the data shown in figure 6.The discharge numbers and DBS system information has been suppressed for clarity.
spectrum resulting from a solution to the inverse problem.The parameters of an analytic δn model are optimized by comparing P inv.s (k meas. ) with the fit to DBS data, P fit s (k) (equation ( 14)) at each k = k meas. .The output of solving the inverse problem is both a set of optimized δn model parameters, and a P inv.
s (k meas. ) spectrum that matches the DBS data (assuming a good fit).We use the following function presented in [11] as our model δn 2 spectrum, This power-law model has five parameters: These parameters are not allowed to vary with radius as this would lead to 'over-fitting'.
Importantly, DBS measurements cannot isolate the effects of the normal (radial) wavenumber spectrum on the scattered power.Therefore, the k n parameters [w n , γ] are underconstrained by the DBS data.Inspired by the TGLF model, we fix γ = 2 to form a Lorentzian in k n .Multiple fixed values of w n were tested.We found the model to be somewhat insensitive to the precise value so long as the spectrum is not highly broadened in the k n -direction, diminishing the cutoff k b resolution.Therefore, we arbitrarily select a scale factor, w n = 2 cm −1 .This simplifies the model allowing us to iterate only the binormal parameters, [k * b , w b , β].For each optimizer iteration, the P inv.
s (k meas. ) spectrum is calculated using a consistent δn.Then, residuals are evaluated (in log-space) according to: Res.= | log P inv.s (k) − log P fit s (k)|.Based on these residuals, the δn 2 parameters ⃗ µ are updated and the process is repeated until convergence.The inverse problem solution, P inv.
s (k meas.), is shown in figure 11.The optimal binormal parameters for the model in equation (15) It should be emphasized that convergence is not unique using this approach, and that the functional form of equation ( 15) injects bias.It was found that convergence is also possible with other models.For example, we also attempted an exponential-type spectrum of the form, Once again, we fix the parameters of the normal-direction such that γ = 1, and w n = 2 cm −1 .The optimal binormal parameters in this case were: Once again, the binormal shift parameter, k * b = 2.4 cm −1 , corresponds to k y ≈ 1.The exponent, β ≈ 1.0 indicates an almost purely exponential δn spectrum.Interestingly, the scale factor, w b , from both the power-law fit and the exponential fit are approximately equal.This is likely driven by the fact that we have fixed w n in both cases (to the same value).The ratio of w n /w b is indicative of the anisotropy in the δn spectrum.

DBS Ps(kmeas.) measurements
The measured DBS P s (k meas. ) spectrum presented in figure 6 exhibits variable spectral decay over the wavenumber range 0.5 ⩽ k ⩽ 16 cm −1 .Table 1 provides a comparison with results from other fusion experiments.Previous studies generally report δn 2 wavenumber spectra with variable decay, including a 'knee' feature near k meas.ρ s = 1, and an intermediate-scale spectral decay ν ≈ 3 (with ν : δn 2 ∼ k −ν ).Our results are similar to what has been previously reported.We observe a spectral 'knee' at k meas.≈ 3.5 cm −1 (k y ≈ 1.4) and an intermediatescale decay index of ν ≈ 2.6.
Variable spectral decay (in log-log space) distinguishes plasma fluctuations from fluid turbulence where a single spectral index may be applicable from driving at a large scale (low k) down to dissipation at a small scale (high k).Plasmas naturally exhibit driving at multiple scales due to an array of unstable modes (e.g.ITG/TEM/ETG).The observed variation in the decay index, ν, could be a symptom of different turbulence regimes, where the final ν value is a nonlinear balance of driving, dissipation mechanisms, and cascade phenomena (such as three wave coupling).
Nonlinear multi-scale gyrokinetic simulations by Howard et al [41] found that the relative importance of ITG/ETG impacted anisotropy in the density spectrum and the spectral energy transfer mechanism.Interestingly, their cases with strong ETG activity exhibited a 'shoulder' in the simulated δn 2 spectrum at k y ≈ 3.5 and a non-local inverse energy cascade from ETG wavenumbers (k y ∼ [2-7.5])extending to ITG wavenumbers (k y ∼ 0.4).
Another notable feature of the P s (k) spectrum shown in figure 6 is the rapid decay (∼ k −9.4 ) at higher-k.Previous studies have also reported spectral decay indices, ν > 7 including several examples from the TJ-II stellarator [35,36].While a value of ν ≈ 9 (or ζ ≈ 4 in the exponential kρ s case) is seemingly large, these are well within the range of previously reported decay indices (exponential spectra are discussed later in this subsection).
DBS mismatch attenuation was also considered when analyzing the seemingly-rapid spectral decay at high k meas. .During the experiment, toroidal steering with the DBS240 system was employed to maintain mismatch-angle θ m < 1 • (at cutoff) for the entire DBS240 dataset.Additionally, the effects of mismatch are included in the SCOTTY calculation of the DBS weighting function.Therefore, any deviations from θ m = 0 are accounted for in the synthetic diagnostic.
The poloidal dependence of turbulence amplitude was also considered in investigating the seemingly-rapid spectral decay at high k meas. .Returning to figure 5(a), it is clear that high-k trajectories tend to probe further away from θ = 0 (θ being the poloidal angle around the flux surface).The highest-k measurements probe closer to θ ≈ 20 • .To investigate whether this could affect the spectral decay, we calculated the theoretical poloidal dependence of the potential spectrum (G(θ) in equation ( 11) [19]) included in the TGLF-SAT2 model.For the plasma considered here, the variation in G(θ) for 0 ⩽ θ ⩽ 20 • was negligible (< 5%).Therefore, within the TGLF-SAT2 theoretical framework, the rapid spectral decay cannot be attributed to the poloidal dependence of turbulence.
In section 6 we made the observation that the measured P s (k) spectrum appears to decay exponentially rather than according to a power-law.Hennequin et al [37] reported collective TS measurements on Tore Supra with an exponentially decaying spectrum, |δn| 2 ∼ e −ζ(k θ ρi) with a value of ζ ≈ 4 over 0.5 ⩽ k θ ρ i ⩽ 2.5 for a variety of plasma conditions.Similarly, DBS measurements by Schmitz et al [42] on DIII-D reported exponentially decaying spectra, |δn| 2 ∼ e −ζ(k θ ρs) , with ζ ≈ 3 in L-mode and ζ ≈ 4 in H-mode for k θ ρ s < 6.
Physically, we can interpret the value of ζ if we take P s ∼ |δn(k)| 2 ∼ S(k) as the Fourier transform of the δn autocorrelation function (cf equation ( 3)).Then, an exponential spectrum in k-space corresponds to a Lorentzian in real-space.If we define the correlation length of δn fluctuations in the  (6,15) <3, 3 Fwd.Collective TS Line-averaged Tore Supra L-mode [40] binormal direction, L b , as the half-width at half-maximum of the Lorentzian autocorrelation function, then L b = ζρ s .Thus, in this simplified picture, ζ is the L b correlation length in units of ρ s .

Synthetic DBS modeling
We have developed a novel synthetic DBS diagnostic to connect the measured P s (k meas. ) with models of the underlying density fluctuation spectrum δn(k).The synthetic diagnostic allows us to approximately separate diagnostic effects from measurements of δn.This is particularly important for DBS given that δn can affect signal localization and the weighting function, Û, can affect the scattered power.We will now discuss further the impact Û has on the relationship between P s and δn.
Figure 10 shows a situation where diagnostic effects are somewhat minimal.The total integrand (figure 10(d)) has roughly the same shape as the probed ⟨|δn| 2 ⟩ along the trajectory (figure 10(c)).Therefore, when P s is calculated by integration, the result is essentially a constant multiplied by the spectrum, i.e.P s ∝ ⟨|δn| 2 ⟩.If this proportionality holds true for every probed wavenumber, then it would be a reasonable approximation to say DBS directly measures ⟨|δn| 2 ⟩ (with some overall normalization).
The analysis presented here does not assume P s ∝ ⟨|δn| 2 ⟩ a priori.Our synthetic DBS diagnostic accounts for variations in the weighting function across the wavenumber scan.To demonstrate instrumental ( Û) effects, we re-calculate P TGLF s (k meas. ) with fixed TGLF-δn from the nominal probed location (ρ = 0.67).This removes the radial variation of the turbulence spectrum from the model.This recalculation is shown in figure 12.In figure 12(a) we show the k n = 0 slice of the TGLF δn 2 spectrum at ρ = 0.67 on the same scale as the P TGLF s .This is essentially a 'before vs. after' view of the convolution with Û(k).Outside of the intermediate-k regime the curve in figure 12(b) deviates from horizontal significantly.This implies that the simple P s ∝ ⟨|δn| 2 ⟩ relationship does not hold across the wavenumber scan.Diagnostic effects hidden inside of the weighting function Û, produce an overall relationship P s ∝ k −0.76 ⟨|δn| 2 ⟩ (as indicated by the linear fit in figure 12(b)).We offer the following explanation: At lowerk meas., the 'ray' term (∼k −2 ) in the P s integrand can begin to out-pace a weakly-decaying δn spectrum.At higher-k meas., where there is little spectral power in the δn spectrum, the ray term once again can begin to dominate the response.Both affects produce an artificially steeper synthetic P s (k).
The differences between the ⟨|δn| 2 ⟩ and P TGLF s in figure 12(a) highlight the importance of treating DBS measurements with a synthetic diagnostic framework.It should be stressed that the precise ratio in figure 12(b) is not a general result.The specifics of the DBS Û vary with the background plasma and the choice of frequency, polarization, and launchposition on the machine.
Finally, note that the P TGLF s (k meas. ) spectrum shown in figure 12(where the radial variation of the δn spectrum is ignored) has worse agreement with the DBS measurements when compared with the P TGLF s (k meas. ) shown in figure 11(where we use cutoff-localized TGLF-δn spectra).In particular, the ρ = 0.67 spectrum deviates at lower-k, and around the 'knee' at k y ≈ 1.This implies that the radial variation in the δn spectrum is a significant factor in the synthetic DBS result.

Scans of TGLF parameters
The results presented in figure 11 show remarkable agreement between the synthetic P TGLF s (k) and the DBS data.While investigating the radial variation in the δn spectra, it was found that these TGLF results are sensitive to the inverse temperature gradient scale-length, R/L Te (TGLF parameter RLTS_1). Figure 13 illustrates the results of an R/L Te -scan carried out by varying RLTS_1 by ±10% for all TGLF radii.
When R/L Te is increased by +10% we observe a substantial rise in the P TGLF s (k meas. ) spectrum over intermediate wavenumbers (figure 13).This change in P TGLF s stems from R/L Te -driven changes in the underlying δn spectra.Figure 14 illustrates the δn 2 spectra at ρ = 0.67 for the nominal R/L Te (figure 14(a)) and +10% R/L Te (figure 14(b)).Figure 14(a) is identical to figure 9 with a log-k b axis.
When R/L Te increases, we observe increased amplitude of both the lower-k (k b ≈ 1.4 cm −1 , k y ≈ 0.5) and the higherk (k b ≈ 2.8 cm −1 , k y ≈ 1.0) modes visible in figure 14(a).Figure 14(b) shows the higher-k mode expands in k-space, becoming the dominant feature in the δn spectrum.Increasing R/L Te also drives a sharp, nonlinear increase in the electron heat-flux, Q e (not shown).The sensitive relationship between Q e and R/L Te suggests transport is driven by TEM/ETG modes.TGLF predicts that given 10% larger R/L Te , the observed P s spectrum would be peaked at k y ≈ 1 (figure 13) as a direct result of changes in the δn spectrum (figure 14).The fact that the DBS measurements do not show a clear peak at k y ≈ 1 suggests that this mode remained subdominant to  lower-k y modes during the experiment or, through a nonlinear process not captured by TGLF, was able to distribute power in k-space.
Meanwhile, when R/L Te is decreased by 10%, figure 13 shows the P TGLF s decreases over intermediate wavenumbers and at the lowest k.When we examine the −10%R/L Te δn(k n , k b ) spectrum (not shown), it appears identical to figure 14(a) with reduced amplitude.In both the nominal and −10%R/L Te cases we observe a possible ITG mode at k b ≈ 1.4 cm −1 (k y ≈ 0.5).The eigenvalue spectrum in figure 8 indicates that this mode moves in the ion diamagnetic drift direction with a growth rate slightly above the local γ E×B shearing rate.

Conclusion
We have presented new measurements of the density fluctuation wavenumber spectrum in the mid-radius (ρ ≈ 0.7) of ECH H-mode discharges using DBS.The back-scattered power wavenumber spectrum, P s (k), was found to have variable spectral decay: roughly flat at low-k (k −0.6 ), moderate decay (k −2.6 ) at intermediate-k, and rapid decay at high-k (k −9.4 ).We also observe a spectral 'knee' at k meas.ρ s ≈ 1, similar to past experiments (see table 1).When the same P s (k) spectrum was plotted on a linear scale against kρ s , it shows variable exponential decay (e −ζ(kρs) ) with ζ ≈ 1.7 for kρ s < 3.5 and ζ ≈ 3.9 for kρ s > 3.5.
We developed a novel synthetic DBS diagnostic combining the beam-DBS model for the weighting function with model δn(k) spectra.The DBS weighting function is calculated using the SCOTTY 3D beam-tracing code.This model includes linear scattering effects important for interpreting DBS signals, e.g.amplification of the probing electric field near cutoff.The SCOTTY code also captures 3D diagnostic effects such as mismatch attenuation.Because the model is linear, it does not include multiple scattering effects nor will it capture the diagnostic response to large amplitude density fluctuations (relative to the background density).
We tested both simulated and analytic models for the δn(k) spectrum using our synthetic DBS diagnostic.In the forward modeling case we used the quasilinear gyro-fluid code TGLF to approximate δn(k).The TGLF model spectrum was then compared with DBS measurements after being convolved with the DBS weighting function.In the inverse modeling case we tested analytic models of the turbulence spectrum by fitting δn(k) model parameters to the DBS measurements.
We found remarkable agreement between DBS P s (k) measurements and the synthetic P TGLF s (k) when TGLF is run at multiple radial locations using the experimental profiles/gradients.The final P TGLF s (k) curve shown in figure 11 includes the radial variation in the turbulence.Scans of physics parameters in TGLF showed that ±10% changes in R/L Te significantly alter the δn spectrum predicted by TGLF.In the +10%R/L Te case, we observe a dominant electron-directed mode at k y ≈ 1 translating to a peak in the predicted P TGLF s .The fact that DBS measurements do not observe this peak, and agree more with the nominal R/L Te case, suggests that these TEM/ETG modes are subdominant in the experiment.
Meanwhile, in the nominal and −10%R/L Te cases we observed a dominant ion-directed mode at k y ≈ 0.5 in the δn spectrum predicted by TGLF.This translated to a maximum in the P TGLF s spectrum at low-k.A holistic combination of the simulation results and DBS measurements supports the conclusion that the core of these ECH H-mode plasma exhibits mixed-mode activity with an ion mode at k y ≈ 0.5 and an R/L Te -driven electron mode at k y ≈ 1 contributing to weak spectral decay at low/intermediate-k.
Future work includes comparing and benchmarking models included in our synthetic DBS diagnostic.Both the DBS weighting function and the δn spectrum can be validated/compared with higher fidelity codes.On the one hand, The TGLF approximation of δn can be compared with nonlinear gyrokinetics simulations.We believe this is an important comparison to make in order to distinguish diagnostic effects from properties of the turbulence, e.g.verifying the TGLF G(θ) model for the poloidal dependence of the turbulence amplitude.On the other hand, the SCOTTY model for the DBS weighting function can be compared with 3D full-wave simulations.While these other methods are undeniably higher fidelity, our synthetic diagnostic presented here has many advantages.Being built on reduced models (TGLF and beam-tracing) substantially reduces development time and enables rapid and informative scans of physics parameters.In addition to exploratory parameter scans, the speed of these methods would allow us to propagate profile uncertainties through TGLF to the final synthetic P s spectrum.
We believe that the work presented here offers a new avenue for testing turbulence models with DBS measurements.If measurements of wavenumber spectra (and the corresponding synthetic diagnostics) are available more routinely, it should be possible to treat the measured spectra as additional constraints for flux-matching workflows such as the TGYRO code [44].Our confidence and predictive capabilities would greatly benefit from simulations simultaneously matching both the transport fluxes and the measured turbulence spectra.
Part of the data analysis was performed using the OMFIT integrated modeling framework [45].Special thanks to J Ruiz-Ruiz for helpful conversations concerning wavenumber metrics.
This material is based upon work supported by the U.S. Department of Energy, Office of Science, Office of Fusion Energy Sciences, using the DIII-D National Fusion Facility, a DOE Office of Science user facility, under Award(s) DE-FC02-04ER54698, and DE-SC0019352.

Disclaimer
This report was prepared as an account of work sponsored by an agency of the United States Government.Neither the United States Government nor any agency thereof, nor any of their employees, makes any warranty, express or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights.Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof.The views and opinions of authors expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof. is adapted from work by Ruiz et al [11].We use the CGYRO code conventions and definitions from Candy et al [20].
The specific relationship between (k x , k y ) and (k n , k b ) in equation ( 10) is formally only applicable to UDS flux surfaces at the low field side midplane (θ = 0).The general equation applicable to arbitrary flux surfaces, and arbitrary θ, originates from the Fourier decomposition of fluctuating quantities in gyrokinetics.For example, the fluctuating density is written as, δn (x, α) = In the UDS and θ = 0 limit, equation (C.3) reduces to equation (10).It's important to note that all DBS measurements presented here come from θ ⩽ 20 • , and from radial locations ρ in the plasma core/mid-radius.To apply equations ( 10) and (C.3) we use the Miller Extended Harmonic parameterization [47] of flux surfaces based on our magnetic equilibrium.
Figure C1 illustrates the approximation made by using equation (10) as opposed to (C.3) at ρ = 0.67. Figure C1(a) compares the magnitude terms in the binormal wavenumber mapping.Specifically, we compare ( b × ∇r/|∇r|) • ∇α with the θ = 0, UDS limit: (1/R 2 + (∂α/∂θ) 2 /(κr) 2 ) 1/2 .Figure C1(b) compares the magnitude of terms in the normal wavenumber mapping.Specifically, we show the θ dependence of |∇r|, and the magnitude of the second term in equation (C.3b) relative to the first term.Generally, we find that using equation (10) for this plasma, over the range of ρ and θ probed by DBS results in an error ⩽5% in converting k x , k y to k n , k b .

Figure 1 .
Figure 1.Various tokamak coordinate systems.(a) illustrates a 3D section of a toroidal flux surface with cylindrical coordinates (R, ϕ, Z) and the simple-toroidal coordinates (r, θ, ϕ).(a) also shows the magnetic field B in purple and a probing DBS ray in red.(b) shows the plane perpendicular to B defining the normal and binormal directions (ên,ê b ).(c) shows the unwrapped flux-surface where B becomes a straight line.
(b)  shows the auxiliary heating P ECH , and P NBI .Figure2(c) shows the core electron and ion temperatures, T e and T i .Figure2(d)shows the D α filterscope.

Figure 2 .
Figure 2. Evolution of DIII-D discharge 189998.(a) shows the line-averaged density (×10 13 cm −3 ) and the plasma current in MA.(b) shows both the ECH and NBI heating in MW.(c) shows the core Te and the core T i in keV.(d) shows the Dα filterscope (a.u.).

Figure 3 .
Figure 3. Radial profiles (a)-(d) of ne, Te, T i , and n Z and their normalized inverse gradient scale-lengths (e)-(h) for discharge 189998, 3000 ms.(a) shows the electron density with data from the Thomson scattering (TS) and profile reflectometry (REFL) diagnostics.(b) show the electron temperature profile with data from TS and electron-cyclotron emission (ECE).The deposition of ECH is also plotted in (b).(c) shows the main-ion temperature measured with main-ion charge exchange recombination spectroscopy (miCER).(d) shows the carbon density, n Z , measured by impurity-CER (impCER).(e)-(h) provide profiles of the inverse gradient scale-lengths normalized with respect to the separatrix minor radius, a.In frames (a)-(d) the solid green lines are fits to the data.In all frames the shaded vertical band indicates the radial region of interest.

Figure 4 .
Figure 4. (a) DBS spectrogram for the 72.5 GHz channel of the DBS240 system for discharge 190129.(b) shows the synchronized Dα filterscope signal.The shaded spans in (b) represent the 50%-95% ELM-phase regions which are considered free from the influence of the preceding ELM and used for further DBS analysis.Several small signals are visible in the Dα trace which were below the threshold used for ELM detection.(c) shows time-averaged spectrograms for three discharges: 190129 (kmeas.≈ 2.5 cm −1 ), 189997 (kmeas.≈ 11 cm −1 ), and 189994 (kmeas.≈ 16 cm −1 ).Negative frequency in (a) and (c) indicates lab-frame motion in the electron diamagnetic drift direction.

Figure 5 .
Figure 5. Ray tracing results illustrating the DBS launch-angle scan.(a) shows a poloidal projection of the 3D ray trajectories for a subset of the scan.The equilibrium plotted in (a) is from discharge 189998, 3000 ms.(b) shows a top-down projection of the ray trajectories in (a).The y-axis of (b) indicates propagation in the toroidal (tor.)direction.(c) plots all probed wavenumbers and radial locations for the entire dataset.The density profile from 189998, 3000 ms is provided against the right-hand axis for reference.(d) shows a histogram of the probed-locations along with a Gaussian approximation with average and standard deviation, ρ ∼ 0.67 ± 0.07.

Figure 6 .
Figure 6.Back-scattered power wavenumber spectrum.The DBS data follows the same marker/color scheme as figure5(c).Black dashed lines represent linear (in log-log) fits with each spectral index provided by annotations.The linear fits, which normally pass through the data points, have been displaced to make them visible.The red line is the result of a nonlinear fit using the function in equation(14).Note the secondary x-axis above the plot providing a mapping between the real-space kmeas., and the normalized ky.

Figure 7 .
Figure 7.Back-scattered power wavenumber spectrum.The DBS data follows the same marker/color scheme as figure5(c).Black lines indicate linear (in semilogy) fits with each exponential decay index provided by annotations.The linear fits, which normally pass through the data points, have been displaced to make them visible.Note, the secondary x-axis above the plot providing a mapping between the experimental kmeas.ρs,and the normalized ky.

Figure 8 .
Figure 8.Comparison of linear eigenvalue spectra calculated by TGLF and CGYRO.Both TGLF and CGYRO are run at ρ = 0.67 using identical inputs from discharge 189998, 3000 ms.The 'x' markers on the CGYRO line indicate cases where long-lived (>300 a/cs), weak oscillations prevented complete convergence.

Figure 9 .
Figure 9. TGLF δn 2 (kn, k b ) spectrum at ρ = 0.67.The red line is a projection of a single DBS trajectory from the wavenumber scan.The arrows indicate the direction DBS traverses through the measured (kn,meas., k b,meas.) phase-space along the path of the beam through the plasma.

Figure 10 .
Figure 10.Components of the beam-DBS model computed by the SCOTTY code.(a) and (b) show the 'Ray' and 'Polarization' parts of the weighting function.(c) shows the δn 2 part of the integrand evaluated along the trajectory of the beam.(d) shows the total integrand (i.e. the product of the preceding terms).The 'mismatch' is not shown as it is roughly a constant between 0.95 and 1.0 for this case.Finally, the shaded region connects this figure with the path shown in figure 9.

Figure 11 .
Figure 11.Comparison of synthetic Ps(k) results with DBS measurements.The forward-model P TGLF s (k) is shown in blue while the inverse-model P inv.s (k) is shown in turquoise.The upper x-axis is the theory ky value (calculated using equation (10a)) while the bottom x-axis is the measured wavenumber in cm −1 .The P TGLF s (k) is labeled 'variable ρ' indicating the radial variation in the turbulence (as predicted by TGLF) is included in the model.The inverse solution is labeled 'equation (15)' to indicate that this P inv.s (k) used equation(15) for the underlying δn spectrum.The DBS data shown in figure11is identical to the data shown in figure6.The discharge numbers and DBS system information has been suppressed for clarity.

Figure 12 (
b) shows the ratio of P TGLF s /⟨|δn| 2 ⟩.If the approximation P s ∝ ⟨|δn| 2 ⟩ were true, the line in figure 12(b) would be horizontal, indicating the constant of proportionality.Interestingly, close inspection of figure 12(b) shows that over the intermediate-k regime (2 < k < 7 cm −1 ) the ratio

Figure 12 .
Figure 12.Illustration of the relationship between Ps and δn.(a) shows the synthetic P TGLF s (k) along with the kn = 0 slice of the underlying δn 2 spectrum (dashed line) and the fit to the DBS measurements (red line).Unlike figure 11, the P TGLF s (k) shown in (a) uses a single radial location (ρ = 0.67).(b) shows the ratio of P TGLF s (k)/δn 2 on a log-log scale along with a linear (power-law) fit.

Figure 13 .
Figure 13.Comparison of P TGLF s (k) spectra for various values of R/L Te .The measurements are represented by their fit (red line).The scan of R/L Te is carried out at all cutoff locations with all other parameters held constant.The overall normalization used in the nominal R/L Te case is also used in the ±10% cases.

Figure 14 .
Figure 14.Comparison of the TGLF δn 2 (kn, k b ) spectra for the nominal case and the case with +10% R/L Te .Both cases use ρ = 0.67.The data shown in (a) is identical to what is shown in figure 9.

Figure C1 .
Figure C1.Comparison of key terms in the general wavenumber mapping (equation (C.3)) with terms in the θ = 0, UDS limit (equation (10)) for the plasma studied in this work (at ρ = 0.67).(a) compares terms in the binormal wavenumber mapping.(b) compares terms in the normal wavenumber mapping.Note that all DBS data comes from θ < 20 • corresponding to θ/π ≲ 0.1 indicated by the yellow shaded region on the left side of both (a) and (b).

Table 1 .
Published wavenumber spectrum results from other fusion experiments.To assemble this table we have equated Ps, δn 2 , and (δn/n) 2 .We have multiplied the spectral-index by two in any cases where authors have published a spectrum ∝ δn.Data shown in the table comes from a variety of confinement regimes, this may have an effect on the range of measured wavenumbers and the decay index, ν.