Scattering of diffracting beams of electron cyclotron waves by random density fluctuations in inhomogeneous plasmas

The physics and first results of the new WKBeam code for electron cyclotron beams in tokamak plasmas are presented. This code is developed on the basis of a kinetic radiative transfer model which is general enough to account for the effects of diffraction and density fluctuations on the beam. Our preliminary numerical results show a significant broadening of the power deposition profile in ITER due to scattering from random density fluctuations at the plasma edge, while such scattering effects are found to be negligible in medium-size tokamaks like ASDEX upgrade.


Introduction
This contribution addresses the physics basis (section 2) and reports preliminary results (section 3) of a new code, WKBeam [1], for electron cyclotron beams in weakly turbulent plasmas.
The quantitative description of the effects of turbulent fluctuations on high-frequency wave beams in tokamaks is not a new problem, and yet it has recently received renewed attention.That was triggered by the work of Tsironis et al. [2] where such effects were found important for ITER on the basis of a simplified but clear physics model.Successive works reported more detailed calculations based upon the application of existing quasi-optical and ray tracing codes [3][4][5], and on ad-hoc methods [6,7].
Contributing to this effort, the WKBeam code relies on the well established framework of kinetic models of radiative transfer [8][9][10] and particularly on the wave kinetic formalism of McDonald [9].The result is a flexible code which accounts for absorption, diffraction and can compute the statistically averaged effect of electron density fluctuations in realistic (experimental or numerical) tokamak equilibria.

The steady-state wave kinetic equation for wave beams
First the theoretical framework for the WKBeam code is outlined with emphasis on physics and basic underlying approximations, while a complete mathematical derivation will be presented elsewhere.

The wave equation
Electron cyclotron beams of angular frequency ω are fully described by the equation for the wave electric field a Corresponding author : Omar.Maj(at)ipp.mpg.deE(ω, x), that is, where the spatial coordinates x = (x 1 , . . ., x d ), d = 3, are normalized to the scale L of the medium inhomogeneity, and κ = k 0 L, with k 0 = ω/c being the wave vector in free space.Here, the response of the medium is accounted for by the linear operator ε, relating the electric displacement εE(ω, x) to the wave electric field E(ω, x).For hot plasmas, such relationship is generally non-local, that is, εE(ω, x) = K(ω, x, x )E(ω, x )dx .
A better characterization of this operator can be obtained by means of a mathematical transformation σ w referred to as (semiclassical) Weyl symbol map [8].The definition of the Weyl symbol resembles a symmetrized Fourier transform of the kernel of the operator, namely, The result is a function of both position x and refractive index vector N.The combined position-refractive-index space (x, N) is the geometrical optics phase space, hence the Weyl symbol map is often regarded as a phase space representation of the operator [8].As an example, one can consider the case of a cold magnetized plasma for which the operator ε reduces to the multiplication by the cold plasma dielectric tensor ε cold (ω, x), so that the kernel is where δ is the Dirac measure which implies that, in the cold limit, the electromagnetic response of a plasma is local.Correspondingly, the Weyl symbol of the operator is readily computed with the result that which shows that, in this case, the Weyl symbol is one and the same with the cold plasma dielectric tensor.Analogously, for a hot (i.e., including temperature effects) uniform plasma, the Weyl symbol coincides with the usual dielectric tensor.For a generic non-uniform hot plasma, a common approximation (the adiabatic approximation) is invoked: The Weyl symbol is identified with the dielectric tensor of a hot uniform plasma ε hot (ω, N) with plasma parameters (magnetic field, electron density, and electron temperature) evaluated at their local value at the position x.
Upon inverting the Fourier transform in equation ( 3) one can write the kernel K in term of the Weyl symbol, and equation (2) becomes which is known as the Weyl quantization of the tensor ε(ω, x, N) and one writes ε = Op w (ε), meaning that the operator ε is obtained from the tensor ε via (4).
For electron cyclotron beams, it is usually believed that the cold plasma approximation is acceptable for the description of the propagation, i.e., for the Hermitian part of the dielectric tensor, while temperature effects must be retained in the calculation of power absorption, i.e., in the anti-Hermitian part of the dielectric tensor.Therefore, it is common practice to simplify the dielectric tensor further, where the cold plasma tensor is used for the Hermitian part (superscript h), and the full hot plasma tensor is retained in the anti-Hermitian part (superscript a) only.The latter is scaled by the inverse of the large beam frequency, i.e., we consider weakly dissipative media.
Under such conditions, approximate solutions of equation (1) with (4) can be conveniently obtained by means of the standard high-frequency (semiclassical) asymptotics, such as the geometrical optics or WKB method (as in ray tracing codes), the paraxial WKB method (as in TOR-BEAM [11]), and the complex eikonal method (as in GRAY [12]).Alternatively, one can describe the propagation, diffraction and absorption of wave beams in terms of a steady-state version of the wave kinetic equation [8].
The validity of the semiclassical approximation is controlled by the small parameter 1/κ.In presence of turbulence, however, short scale variation of the plasma profiles can occur and the inhomogeneity scale length L, estimated from fluctuating profiles, can be significantly reduced, thus spoiling the semiclassical approximation.
Nonetheless, if one is interested in statistically averaged effects, it was proven by McDonald [9] that semiclassical methods can still be applied by adapting an argument of Karal and Keller [13] to the wave kinetic equation.
According to such an approach, the background equilibrium is separated from the fluctuations and the scale L is always computed from the background so that κ has its usual large value.

Statistical description of turbulence
In turbulent plasmas, the electron density, temperature and magnetic field are fluctuating quantities.In this work, we study the effect of density fluctuations only, since that is usually considered the most detrimental for tokamak applications [2].Temperature fluctuations cannot contribute to the lowest order equations as they can affect only the anti-Hermitian part of the dielectric tensor (5), which is already a first order correction.Magnetic field perturbations are usually weak, although Balakin and Gospodchikov [4] claimed that, in principle, they could have an effect by changing the group velocity even without significant change of the wave vector.
A first observation is that, strictly speaking, a turbulent plasma is not a stationary medium and the use of the frequency-domain equation (1) may appear inappropriate.
With reference to heating and current drive scenarios, however, the beam is maintained in the plasma for a time interval much longer than the correlation time of the turbulence, and we are interested in the time-average effect of fluctuations on the beam.We assume (and do not discuss this assumption in details) that such time average can be replaced by an ensemble average over stochastic perturbations of plasma parameters, particularly the electron density, in equation (1).
Mathematically, the total electron density is modeled as the sum of the nominal equilibrium profile n e (x) plus a stochastic perturbation where μ is a random field sampled from a probability measure.Here, the strength of the perturbation is assumed to be of order 1/ √ κ; this apparently artificial assumption is common in radiative transfer modeling [10], but can be justified on the basis of mathematical arguments only (the frequency of the beam, κ, and the strength of the fluctuations, δn e , are physically independent parameters).We also assume that the average effect of the fluctuations does not change the equilibrium, namely, where E is the expectation value operator (average) of the random field.A possible non-zero average could be dealt with by redefining the equilibrium.

EPJ Web of Conferences
Correspondingly, the perturbation of the Hermitian part of the dielectric tensor (5) takes the form where I is the identity tensor.This form of the perturbation stems from the fact that the cold dielectric tensor is linear in the density.The perturbation of the hot anti-Hermitian part of ( 5) is ignored for simplicity, because, in view of (6), it amounts to an O(κ −3/2 ) correction to the wave equation and that does not affect the lowest relevant orders in the semiclassical limit κ → +∞.Hence, the dominant perturbation of the dielectric tensor is independent of N, and this simplifies greatly the analysis [9].Equation ( 7) implies E(δε) = 0. Physically, this framework can be interpreted as follows.Let us imagine we know exactly the electric field E(t, ω, x) of a beam in the presence of time-dependent density fluctuations δn e (t, x).In order to compute the timeaverage, we can extract a sample E j (ω, x) = E(t j , ω, x) for each observation time t = t j and compute the average of the sampled time slices.When the typical time scale of turbulent structures is much larger than the propagation time of the beam (given by the size of the domain divided by the group velocity), we can assume that the density fluctuations are frozen to the value δn e (t j , x) (this is known as frozen turbulence approximation) and the electric field E j (ω, x) is well approximated by a solution of the frequency domain equation ( 1) perturbed by δn e (t j , x).If we can exclude any correlation between the turbulence screen-shots at different time, i.e., if Δt = t j+1 − t j τ corr with τ corr the correlation time of the turbulence, each δn e (t j , x) can be regarded as an independent realization of the random field ( 6) and E j (ω, x) can be regarded as the corresponding sample of solutions of equation ( 1) with the stochastic perturbation (8).Then averaging over such time slices t j should be equivalent to the ensemble average over many independent realizations of the random field μ.This heuristic argument also suggests that the scattering of a beam from a coherent structure might not be correctly described in this statistical framework as different time slices would not be statistically independent.A better and more detailed physics analysis of the statistical approach has been given by Peysson and co-workers [5].
Since, in turbulent plasmas, the wave equation ( 1) exhibits a stochastic perturbation, the solution E(ω, x) itself is a random field and we consider the two point correlation matrix, where EE * on the left-hand side is a direct product (dyad), so that C is a tensor.If E is locally squared integrable, the trace of C is proportional to the averaged energy carried by the wave electric field, namely, where Ω is an arbitrary spatial domain and TrC is the trace of the matrix C. On the lines of the approach of McDonald [9], which in turn, is adapted from Karal and Keller [13], we shall obtain an equation for C by means of a careful average combined with the standard Born approximation.The validity of the Born approximation is therefore crucial for the application of this framework.

Wigner matrix and Wigner function
Actually, rather than in the correlation function C(ω, x, x ) itself, we are interested in its phase space representation.This can be readily obtained by the application of the Weyl symbol map, cf.equation ( 3) with K replaced by C, where the two-point correlation function C(ω, x, x ) is regarded as the kernel of the integral operator Ĉ called correlation operator [8].
The Weyl symbol W(ω, x, N) of the correlation operator is referred to as the (averaged) Wigner matrix of the field E(ω, x) and it is the central quantity in our study.
The Wigner matrix and the correlation matrix, being related by a Fourier transform, encode the same information; thus, the Wigner matrix can be used to compute physical observables that depend quadratically on the field.E.g., the energy carried by the electric field in a domain Ω is proportional to TrW(ω, x, N)dxdN.
Let us introduce the notation where |Ω| is the volume of the domain Ω and a(ω, x, N) is a generic matrix-valued function over the phase space.
When Ω = Ω x is chosen to be a small neighborhood of a point x the latter result can be used to compute the spatial and ensemble average of the electric field energy density, The power deposited in a volume Ω can be obtained in the same form; the Poynting theorem for equation (1), namely, shows that the average power deposited in a volume Ω is proportional to Tr(ε a hot W)dxdN. Here, where the numerical factor is obtained upon restoring the physical dimensions and ε a hot /κ is the full (i.e., non rescaled) anti-Hermitian part of the dielectric tensor, cf.equation (5).
In addition to the possibility of computing relevant physical observables, the Wigner matrix has its own physical interpretation.In order to make that clear, let us consider, in one spatial dimension, a scalar field of the form u(ω, x) = e −x 2 /2 2 cos(κN 0 x), (13) which consists of two counter-propagating Gaussian packets with refractive index vectors N = N 0 and N = −N 0 , respectively.The corresponding scalar Wigner function w(ω, x, N) can be computed analytically from equation (10) and it is shown in figure 2. In the spatial variable x, w(ω, x, N) follows the Gaussian envelope of the squared field u 2 ∼ e −x 2 , while in the wave-vector variable N one can clearly distinguish the two spectral contributions of ±N 0 as well as the interference patter at N = N 0 + (−N 0 ) = 0. Away from the interference region where it can attain negative values, the Wigner function w behaves as an energy density in phase space, hence, it is called quasidensity.
As κ → +∞, the oscillations in the interference region increase in frequency but maintain the same amplitude, therefore the high-frequency limit of the Wigner function we are interested in does not exist in strong sense.Nonetheless, regarding w as a measure means that it makes sense only when integrated over a phase space volume as used in equations ( 11) and ( 12), thus, we do not actually need the strong limit.In addition, the stationary phase lemma ensures that the interference pattern will not contribute to the integral of w, so that the high-frequency (semiclassical) limit exists in this weak sense and gives a positive definite energy density over the phase space.Such an energy density is the physical quantity computed by the WKBeam code.

The approach of Karal and Keller and the Born approximation
The wave equation (1) along with (4), (5), and (8) can be written in the compact form where D = Op w (D), is the standard plasma dispersion tensor, and δε is the stochastic perturbation defined in (8).The Born approximation, as it is used here, consists in treating δε as a perturbation (consistently with the ordering ( 6)).Therefore, the wave electric field is written as a formal asymptotic series and the sequence E n is determined by the hierarchy The lowest order term E 0 is a deterministic field corresponding to the solution in absence of fluctuations, while E n are random fields due to the presence of the stochastic perturbation δε.

EPJ Web of Conferences
It is enough to solve the equations for the fluctuating field corrections E n , n ≥ 1, to the lowest order in the semiclassical limit κ → +∞.With this aim, it is possible to use symbol calculus and construct an operator Q such that that is, Q is an approximate inverse of D. Such operator Q is constructed analytically from the Hermitian part of the plasma dispersion tensor (14).This is almost as good as an exact inverse: If κ is large enough, (I + κ −1 R) −1 exists and can be constructed by the Neumann series n (−1/κ) n Rn .Then, the exact recursion relation In particular, we have where the ordering (6) has been accounted for; one also notices that ) and E(E 1 ) = 0.
Here semiclassical calculus is applied to the operator D which is determined by the plasma background and does not include the fluctuations.The short-wave-length limit therefore applies.
The two-point correlation matrix amounts to ), and upon accounting for (16) and multiplying on the left by D, we obtain the equation for the correlation operator in the form When mapped into the phase space by the application of the Weyl symbol map σ w , this equation for the correlation operator yields a corresponding equation for the Wigner matrix, namely, where the Poisson brackets and g are extended to matrices component-wise, namely, {A, B} i j = k {A ik , B k j }, and S(Γ, W) is an integral operator acting on the Wigner matrix W and depending on which is the Wigner function of the density correlation.Particularly, turbulence enters equation (18) through Γ, which is a statistical property and, as such, is expected to exhibit spatial variations on the same scale as equilibrium profiles.Hence, formal semiclassical calculus (highfrequency limit) applies to (18).

The wave kinetic equation for coherent beams
The formal asymptotic analysis of equation (18) in the semiclassical limit κ → +∞ shows that where w α are the two scalar Wigner functions for the ordinary mode (α = O) and the extra-ordinary mode (α = X), and e α (ω, x, N) is the corresponding polarization unit vector in the cold plasma approximation.Hence, to the lowest relevant order in the semiclassical limit, the Wigner matrix W is diagonal on the same basis as the Hermitian part of the plasma dispersion tensor ( 14) and its eigenvalues are the Wigner functions of the two modes supported by the plasma [8].
For the Wigner functions w α , α = O, X, one finds the constrained phase space transport equation S αβ (Γ, w α , w β ), which is the steady-state form of the wave kinetic equation.Here, for each mode α, H α = e * α D h e α is the geometrical optics Hamiltonian, which depends only on the Hermitian part D h of the plasma dispersion tensor ( 14) and describes wave propagation, γ α = e * α ε a hot e α is the absorption coefficient, which dissipates the electric field energy density, and S αβ is the scattering operator, which depends on the Wigner function (19) of the density correlation and describes the average energy exchange with the mode β due to the effects of turbulence.In addition, accounting for the spectral width of the beam, the wave kinetic equation properly describes diffraction effects [14], in spite of being based upon the same Hamiltonian structure as the standard ray tracing.
On one hand, the derivation of (20) follows from the direct application of McDonald theory [9].On the other hand, one should notice a few differences that are peculiar to the present steady-state form.First, the presence of the constraint H α w α = 0, which expresses the fact that there is no wave energy density, i.e., w α = 0, outside the dispersion surface H α = 0.As a consequence, w α is a singular distribution: w α ∝ δ(H α ).Moreover, the Hamiltonian H α is the eigenvalue of the Hermitian part of the dispersion tensor (14), with the frequency ω of the beam being a parameter, and the transported quantity is the Wigner function which is related directly to the electric field energy density.This should be compared with the standard timedependent form of the kinetic equation, where the dispersion relation is solved for the frequency, which plays the role of the Hamiltonian, and the wave action density is the transported quantity.Correspondingly the absorption coefficient γ α as defined here accounts for the dissipation of the electric field energy density rather than of the wave action density.
For the solution of the kinetic equation ( 20) is crucial to have the proper boundary conditions: One should be careful to retain the finite spectral width of the beam in the boundary data, as that contributes essentially to the diffraction of the beam.The consistent procedure to build boundary data for (20) is inferred from the definition of the Wigner function: At the launcher we assume that the wave electric field is given in the form of a Gaussian beam with prescribed beam parameters; the Wigner function of the initial Gaussian field is computed analyt-EC18 01002-p.5 ically yielding the boundary datum w α|Σ (ω, y, N y ) on the plane Σ of the launching mirror, which is parametrized by two-dimensional coordinates y = (y 1 , y 2 ), with N y = (N y,1 , N y,2 ) being the corresponding components of the refractive index vector.Hence, w α|Σ is defined in a fourdimensional hyperplane in the phase space, and, together with the local dispersion relation, gives the boundary data.

The scattering operator
The derivation outlined above gives the scattering operator which expresses the (normalized) rate-of-change of the energy w α (ω, x, N)dxdN contained in an infinitesimal phase space volume centered on (x, N) as a result of cumulative energy gain via scattering N → N, from any N , (first term on the right-hand side) and cumulative energy loss via scattering N → N , to any N , (second term on the right-hand side).It is worth noting that the obtained scattering process involve solely the refractive index vector, i.e., the process is local in space.The differential scattering amplitude σ αβ is found to be where κΓ = O(1) for κ → +∞ in view of (6).Apart from the factor weighting the scattering amplitude with the effect of polarization, we notice that a scattering process N → N is possible only if the target N satisfies the dispersion relation (cf. the factor δ(H α )) and the spectrum of turbulence does not vanish at ΔN = N − N (cf. the factor Γ(ω, x, N − N )).The latter condition means that the scattering process is a momentum conserving wave-wave interaction, namely, where the variation of momentum ΔN is carried by the turbulent fluctuations that act as a momentum reservoir.
(Here we think of the refractive index vector as the momentum in the sense of geometrical optics.)The projection of this identity onto the direction parallel to the magnetic field gives as in tokamak plasmas the spectrum of the turbulence is such that ΔN ≈ 0 since turbulent structures are elongated along the magnetic field lines.Thus, the parallel refractive index does not change much in the scattering process.
In addition, both N and N have to satisfy the dispersion relation of the corresponding mode, which in a magnetized plasma has the form Figure 3. Energy distribution of the beam in the angle φ N for the case of the full scattering operator (green curve labeled "transport") and its diffusion limit (red curve labeled "diffusion"); for reference the evolution without scattering is also shown (blue curve).For parameters in the range of ASDEX upgrade, the diffusion limit overestimates scattering effects (a), while it is a good approximation for large devices like ITER (b).
where N ⊥ = |(I − bb)N| is the refractive index perpendicular to the local magnetic field unit vector b, and If we consider scattering of radiation in the same mode, the dispersion relation fixes N ⊥ , and we can conclude that in magnetized plasmas the scattering process is essentially a change in the orientation angle φ N defined by N = (N ⊥ cos φ N , N ⊥ sin φ N , N ) in an orthogonal frame where b is directed along the third axis.
It is instructive to compute the diffusion limit of the scattering operator (21) for the process O → O and for perpendicular propagation.As a model for Γ, we set Γ ∝ δ(N − N ) in the parallel direction, so that scattering preserves N exactly, and consider a normal distribution of amplitude δn 2 e,0 and spectral width 1/k 0 L ⊥ in the perpendicular direction; these is equivalent to the conditions considered by Tsironis et al. [2], and we expect to obtain the same result.Indeed, the scattering operator for the Omode is approximated by a diffusion in the angle φ N and we find the effective diffusion coefficient [1], If L ⊥ is identified with the blob size, equation ( 25) is the same as the result of Tsironis et al. [2], apart from the numerical factor, which depends on the details of the model, that is, either Γ is the present statistical approach or the shape of Gaussian blobs in reference [2].This proves that the scattering operator obtained by Weyl calculus captures the expected physics in this simple case.Figure 3, however, shows that the diffusion limit is not always valid and, thus, the full scattering operator should be employed.The full scattering operator (21) accounts for crosspolarization scattering as well, i.e., one can have transfer of energy between different wave modes, namely, (β, N ) → (α, N).Again N and N must satisfy the dispersion relation (24) of the corresponding mode, while momentum conservation (23) and the triangle inequality give Since N ≈ N we see that, for cross-polarization scattering to be significant, the distance between the dispersion EPJ Web of Conferences 01002-p.6 curves of the two modes must not exceed the average ΔN of the turbulence.Figure 4 shows the scattering amplitudes integrated in N for the O → O and the O → X scattering processes in ITER, computed along a geometrical optics ray (in the O-mode) and plotted as functions of the normalized poloidal flux ρ, for the case in which the density fluctuations are centered in a narrow layer around the separatrix where the modes are separated.One can see that, in this case, cross-polarization scattering is suppressed since the dispersion curves are sufficiently separated one from the other.
In the WKBeam code, cross-polarization scattering is neglected; however, we retain the possibility to compute the integrated scattering amplitude (as in figure 4) in order to check the validity of this assumption.

Numerical results and conclusions
WKBeam is a Monte Carlo code designed to solve the steady-state form (20) of the wave kinetic equation, neglecting cross-polarization scattering and using the boundary conditions described in section 2.5.The algorithm developed for WKBeam has been rigorously justified but we shall give here a simple physics description only.
The main computation is very similar to ray tracing.Rays are traced and random scattering events are generated along each ray according to the diagonal elements σ αα of the scattering amplitude (22).If τ is the parameter along each ray, then the interval Δτ between two scattering events is drawn from an exponential distribution.Each scattering event changes the phase space position of a ray according to (x, N ) → (x, N), where the "scattering kick" N → N follows the physics outlined in section 2.6.The result is a set of ray trajectories that sample a particular Poisson process with the property that its probability distribution is related to w α .Physical observables of the beam can then be evaluated by a suitable Monte Carlo estimator, e.g., of integrals (11) for the electric field energy density and (12) for the power deposition profile.
The code differs from a standard ray tracing with scattering in the initialization of rays in the phase space.As mentioned in section 2.5, boundary conditions should be constructed from the Wigner function of the launched electric field and that is defined on a four-dimensional hyperplane.The initial ray positions should sample the Wigner function in this four-dimensional space.The number of required rays is therefore greatly larger that a standard ray tracing calculation (typically ∼ 10 5 rays for a converged calculation of power deposition profiles in ITER).This is the computational price we pay for a consistent description of diffraction and scattering.
A detailed verification of the code has been carried out in absence of turbulence: The analytical solutions of all the standard test models for propagation (free space, lens-like medium, and linear layer medium) have been successfully reproduced.
In tokamak geometry and without fluctuations, the power deposition profile for an ASDEX upgrade shot has been compared with the corresponding calculation from the paraxial WKB code TORBEAM [11] which can account for the full tokamak geometry as well as diffraction effects.The agreement in figure 5 suggests that both geometry and diffraction effects are properly dealt with in WKBeam.Verification and benchmarking efforts in presence of fluctuations are in progress; this is more difficult as analytic solutions are not easily available.A more subtle theoretical issue that should be addressed is the assessment of the validity conditions of the Born approximation [15], cf.equation (15), which is a fundamental step in the derivation of the wave kinetic equation.
Even though the verification and benchmarking of scattering processes in WKBeam is not yet complete, first results on the power deposition profiles show an interesting physical difference between present day medium-size tokamaks such as ASDEX upgrade and large machines, ITER in particular.
Figure 6 shows power deposition profiles for the AS-DEX upgrade shot #25484 obtained by WKBeam for different thickness of the turbulence layer at the plasma edge.In all cases fluctuations produce long tails of the profile without affecting the bulk: This modest effect is consistent   with the transport regime predicted for medium size tokamaks, cf.figure 3a, and with the empirical observation of power deposition in present day machines.
For the case of ITER shown in figure 7, the regime is found to be diffusive, cf.figure 3b, with much more important consequences.One should also note that in this case the dependence on the parameters of the turbulence model, i.e., of the Wigner function Γ in equation (19), is much stronger: The development of an accurate and vali-dated model for such quantity appears therefore crucial for reliable predictions.The consequence of these results for the ITER upper launcher have been discussed by Poli et al. in these proceedings.
In summary, the new wave propagation code WKBeam has been developed on the basis of a well-established theoretical framework which relies on the Born approximation and semiclassical limit.The main feature of the code is the possibility to describe, within the limits of such approximations, the average effects of density fluctuations (including short scale fluctuations) together with diffraction effects with either experimental or numerical tokamak equilibria.Preliminary results show that turbulence has a modest effect for medium size devices (as expected), but may be important and deserve in-depth consideration for large devices such as ITER.
is the normalized Poynting flux for harmonic fields and equation (5) has been accounted for.By selecting a volume shell Ω = Ω ρ around EC18 01002-p.3

Figure 1 .
Figure 1.Spatial domain Ω ρ used to compute the power deposition profile (12) in tokamak geometry for the case of an ITER equilibrium.The domain is centered around a flux surface (red) and bounded by two nearby flux surfaces (blue) defining a shell of thickness Δρ in normalized poloidal flux.

Figure 2 .
Figure 2. Contours of the Wigner function w(ω, x, N) in the twodimensional (x, N) phase space for the case of two interfering Gaussian wave packets (13), with N 0 = 0.5 and κ = 10.

Figure 4 .
Figure 4. Total scattering amplitudes for the O → O (blue curve labeled Σ OO ) and for the O → X (green curve labeled Σ OX ) processes.The latter has been multiplied by 10 3 to make it visible on the same scale.For reference the dispersion curves of the Oand X-mode are also shown.The shadowed area represents the average |ΔN| of the turbulence.This is for an ITER equilibrium.

Figure 5 .
Figure 5. Beam and power deposition profile for the shot #25485 in ASDEX upgrade.The power deposition profile, in particular, is in good agreement with the corresponding TORBEAM calculation.The error bars in the WKBeam curve represent the variance of the Monte Carlo estimator for integral (12).

Figure 6 .
Figure 6.Power deposition profile for ASDEX upgrade shot #25485 computed by WKBeam in presence of edge density fluctuations.A scan in the thickness Δρ F in normalized poloidal flux ρ of the turbulence layer around the separatrix is shown together with a reference calculation without scattering effects.One can see that the profiles are consistent with the transport regime of figure 3a and thus are not severely affected by fluctuations.Here, Γ is an anisotropic normal distribution in ΔN = N − N , localized in ρ around the separatrix with width Δρ F , with maximum such that δn e /n e = 0.1, large parallel correlation length L , and perpendicular correlation length L ⊥ = 0.4cm.

Figure 7 .
Figure 7. Power deposition profile for ITER computed by WK-Beam for a scan of turbulence parameters as in figure 6.For ITER the diffusion regime of figure 3b sets in and the effects of scattering appears much more relevant.The same model as in figure 6 is used for Γ, and again δn e /n e = 0.1, but the perpendicular correlation length is L ⊥ = 2cm.