1 Introduction
Radio frequency (RF) waves are extensively used for plasma heating, current drive and suppression of neo-classical tearing modes. RF waves have to propagate through frequently turbulent plasma edge before they reach their target. The edge is populated by filamentary and blob-like structures of different shapes and density contrasts (see Krasheninnikov Reference Krasheninnikov2001; Grulke et al.
Reference Grulke, Terry, LaBombard and Zweben2006; Myra et al.
Reference Myra, D’Ippolito, Stotler, Zweben, Leblanc, Menard, Maqueda and Boedo2006a
; Myra, Russell & D’Ippolito Reference Myra, Russell and D’Ippolito2006b
; Zweben et al.
Reference Zweben, Boedo, Grulke, Hidalgo, Labombard, Maqueda, Scarin and Terry2007; Pigarov, Krasheninnikov & Rognlien Reference Pigarov, Krasheninnikov and Rognlien2012). As a result, the plasma permittivity of these structures differs from the permittivity of the background plasma. This leads to a modification of the propagation characteristics of the incident RF beam. Thus, an estimate of an effective permittivity that can successfully characterize the turbulent edge as a whole is necessary. Furthermore, due to the statistical nature of the turbulent structures (shape, size, contrast), the estimation method must be such as to be easily amenable to statistical treatment. Briefly, the overall endeavour consists of the following three steps. The first step is the subject matter of the present work. In the second step the probability distribution of the geometrical characteristics, the contrast and the respective filling ratios are suitably provided (as close to experimental observations as possible) and a statistically averaged homogenized permittivity is estimated that is as general as possible. In the last step the latter is used for a turbulent layer populated by blobs. In a scheme which involves the so-called
$4\times 4$
matrix technique (Berreman Reference Berreman1972) for a beam, the beam can easily be constructed as a superposition of as many planar waves as necessary. An alternative method for the propagation of a beam through a turbulent environment has also been developed by Snicker et al. (Reference Snicker, Guidi, Kohn, Maj, Weber and Poli2016).
In the present work we develop a homogenization method capable of encompassing the geometrical characteristics of an aggregate of identical structures (as far as their shape, size and contrast are concerned) immersed in the ambient magnetized plasma of the edge. The blobs are ellipsoidal structures of given size of their three independent semi-axes and density contrast, aligned along the ambient magnetic field and immersed in the ambient magnetized plasma that is considered cold. The filling ratio at this stage is a free parameter. The outline of the present paper will be as follows. In the first section, the generalized dyadic Maxwell equation will be presented, along with the explicit laws governing plasma medium, blobs and, finally, the composite medium. In the second section, reference to the previous homogenization approaches will take place, showcasing two main methods (one older and one much more recent). The third section refers to the mathematical formalism of the proposed method as well as results for electron cyclotron plasma and comparison with preexisting ones. In the last section we summarize our results and briefly discuss future extensions.
2 The formalism of homogenization methods
The problem of a plasma medium with embedded blobs can be described as follows. Considering an electron–proton plasma, the protons propagate at much lower velocities than the electrons, due to their higher mass. This will lead them to gather forming filamentary structures. These structures will be elongated toward the external magnetic field. These inclusions are called blobs and, in general, their shape will be ellipsoidal.
Following chap. 2.3 of Mackay & Lakhtakia (Reference Mackay and Lakhtakia2015), after a temporal Fourier transform of the Maxwell equations, the correlation of the electric displacement and magnetic field to the electric field and magnetizing field can be given by:
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn1.gif?pub-status=live)
Where
$\unicode[STIX]{x1D649}_{6\times 6},\unicode[STIX]{x1D647}$
are
$6\times 6$
tensors, with the former dependent on spatial partial derivatives and the latter constant, computed by the frequency and external magnetic field of the plasma medium. Finally,
$\boldsymbol{f}$
is a 6-vector dependent on position, comprising of normalized electric and magnetic field components.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn2.gif?pub-status=live)
Every quantity used in the present notation has been appropriately scaled to be dimensionless. This
$6\times 6$
dyadic notation can be further simplified in the case of a dielectric plasma medium. This is because the magnetic field components remain unaffected in this equation. Thus, the equation which will be of importance here will be
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn3.gif?pub-status=live)
Where now
$\unicode[STIX]{x1D649}$
is the part of the
$6\times 6$
tensor that acts on the electric field, and the plasma dielectric tensor
$\unicode[STIX]{x1D646}$
is given by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn4.gif?pub-status=live)
Definition of another auxiliary set of functions needs to take place to abbreviate calculations. The medium in question is a dielectric plasma with embedded blobs. Let
$\unicode[STIX]{x1D646}_{P},\unicode[STIX]{x1D646}_{B},\unicode[STIX]{x1D646}_{H}$
denote the dielectric tensors of the plasma, blob and composite medium respectively. The composite medium is the actual medium, which includes both the plasma and the ellipsoidal blobs. The aim of the present paper is, of course, the computation of the composite tensor. The volume fraction of the blobs inside the plasma is denoted by
$\unicode[STIX]{x1D70E}$
. Next, it is proper to define the electric fields
$[\boldsymbol{e}]_{P},[\boldsymbol{e}]_{B}$
as the solutions of the differential equations
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn5.gif?pub-status=live)
Meaning they are the equivalent electric fields, had the medium been homogeneous, comprising only of the plasma material in the first case, and blob on the second. Equation (2.3) is transformed under this notation as follows
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn6.gif?pub-status=live)
Define the Green’s functions for plasma and blob as solutions to the equations
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn7.gif?pub-status=live)
where
$\unicode[STIX]{x1D6FF}(\unicode[STIX]{x1D746}-\unicode[STIX]{x1D746}^{\prime })$
is the three-dimensional Dirac delta function, and
$\unicode[STIX]{x1D644}$
is the 3 by 3 identity matrix. The equation that completed the set of equations to be used in previous homogenization formalisms is:
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn8.gif?pub-status=live)
Where
$V_{E}$
is the volume of the ellipsoidal inclusion. Using (2.5) and (2.8) is the basis for previous methods, as will be presented in the next subsection. Note that compared to the present method, previous ones made no usage of the fact that the tensor
$\unicode[STIX]{x1D646}_{H}$
satisfies (2.6).
2.1 Previous formalisms
2.1.1 Depolarization dyadic
Widely used in the context of homogenization formalisms is the notion of the depolarization dyadic, which makes its introduction a priority before continuing. The depolarization dyadic is defined by the following relation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn9.gif?pub-status=live)
In general, the Green’s function of the plasma is unknown, but its Fourier transform is readily extracted by Fourier transforming (2.5). Depolarization dyadics are only computed via this Fourier transform. The first approximation used is the Rayleigh approximation: assuming uniform distribution of the charge density
$\unicode[STIX]{x1D64C}_{\text{equiv}}(\unicode[STIX]{x1D746})$
throughout the blob, equation (2.12) becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn10.gif?pub-status=live)
This approximation leads to the long wavelength approximation: the blob is supposed to be sufficiently smaller than the incoming beam wavelength in order for the Rayleigh approximation to be satisfied.
Utilizing Fourier transforms of Green’s function and of the source term in (2.10) the resulting equation for the depolarization dyadic is (see Mackay & Lakhtakia Reference Mackay and Lakhtakia2015)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn11.gif?pub-status=live)
Where
$q$
is the wavenumber and
$\unicode[STIX]{x1D650}$
is a dyadic that transforms spherical coordinates into elliptical ones. The dyadic itself can be approximated, because its numerical computation is cumbersome. It can be separated into two different dyadics, the one utilizing the Green’s function limit as the wavenumber
$q\rightarrow \infty$
and the second one being the full dyadic minus the first approximation. There is a reason for such a method: the full dyadic may initially produce many errors when computed, due to same order of
$q$
-polynomial in the numerator and denominator of the Green’s function. However, these two separate dyadics are smoother to extract numerically, because the first Green’s function is
$q$
independent, and the second one has a numerator of smaller
$q$
order than the denominator.
To compute the dyadic more efficiently, these equations are being used:
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn12.gif?pub-status=live)
2.1.2 A simple mixing formula
A very simple mixing formula for a dielectric material with embedded blobs was the method described at Sihvola (Reference Sihvola1996). Using the simplest depolarization dyadic approximation, the resulting equation is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn13.gif?pub-status=live)
This mixing formula can of course be improved by using the full depolarization dyadic. An important aspect of the above equation is its symmetry: no inversion of parameter definition needs to take place when the volume fraction exceeds 50 % because at volume fraction 0 and 1 the results are as expected (plasma tensor and blob tensor respectively).
2.2 Differential Maxwell–Garnett formalism
Homogenization can also happen incrementally, in the following sense. Homogenizing a material while only incorporating one blob at a time, finding a new composite material for each step. It is obvious that the number of steps needed to reach the final material is too large, enough to be considered
$N\rightarrow \infty$
. Thus, it can be replaced with the differential Maxwell–Garnett formalism (Mackay & Lakhtakia Reference Mackay and Lakhtakia2015). The differential equation to be solved is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn14.gif?pub-status=live)
This can be approximated numerically. Results depend on the choice of depolarization dyadic, as well as wavelength of the incoming beam and volume fraction of the blobs. In cases where the equation cannot be solved analytically, it is essentially the same as the incremental method in terms of computational efficiency, while retaining limitations of previous formalisms.
3 The proposed new method
The limitations of methods illustrated above are the long wavelength approximation and behaviour of solutions for large filling ratios. The Rayleigh approximation assumes that the transverse blob dimensions are small compared to the wavelength of the incoming electromagnetic wave. In the scrape-off layer this approximation is difficult to satisfy. For example, for electron cyclotron waves, the wavelength is of the order of millimetres while the blob radius is of the order of a few centimetres. The Rayleigh approximation stems from the fact that, usually, regarding the depolarization dyadic, only the first dyadic in (2.12) is used in current formalisms (trading, of course, physical information for better computational efficiency). However, only the second one incorporates the size of the blob, as well as its shape. Our own method makes no usage of the long wavelength approximation. To illustrate the method that is to be used in the present paper, two new quantities need to be introduced. For the resulting equation to be dimensionless, define the normalized wavenumber
$q=kc/\unicode[STIX]{x1D714}$
, where
$k,c,\unicode[STIX]{x1D714}$
are normal wavenumber, speed of light and angular frequency respectively. Next, the dyadic that is to be used when changing coordinate variables from spherical to ellipsoidal is denoted by
$\unicode[STIX]{x1D640}$
which consists of the ellipsoid semiaxes
$a,b,c$
(their dimensions are length dimensions, however the dyadic is normalized to be dimension free).
Fourier transforming equation (2.5) and (2.6) gives
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn15.gif?pub-status=live)
where the matrix
$\unicode[STIX]{x1D648}=\boldsymbol{q}\boldsymbol{q}-q^{2}\unicode[STIX]{x1D644}$
. Equations (3.1) can be solved to provide the normalized Fourier transformed electric fields as functions of
$\boldsymbol{q}$
and the dielectric tensor components. However, the homogeneous nature of the above equations leaves a ‘freedom of choice’. That is, one can multiply the solutions by any scalar function and maintain validity of the equations. This freedom can be used to ensure that boundary conditions for the electric fields hold. However, even the simple inverse square law as distance goes to infinity can be enough to dictate choice for that function, and that is going to be used in the present work.
3.1 Homogenization equation
The starting point is that an exact analytical solution for (2.8) is needed, at first for the case of a single blob to provide the basis upon which our method will continue. The general solution is given in the following form:
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn16.gif?pub-status=live)
The proof is straightforward (Mackay & Lakhtakia Reference Mackay and Lakhtakia2015) and it can be facilitated by using (2.5) and (2.7).
It must be emphasized that (3.2) does not depend on any approximations so far. However, if one approximates
$\unicode[STIX]{x1D64C}_{\text{equiv}}(\unicode[STIX]{x1D746}^{\prime })$
by
$\unicode[STIX]{x1D64C}_{\text{equiv}}(\unicode[STIX]{x1D746})$
instead, then this is equivalent to assuming that there is no equivalent charge variation over the wavelength of the wave. This is, of course, the Rayleigh approximation and is the basis of all previous studies. We will not be making this approximation and consequently we will not impose any restrictions on the ratio of the wavelength to the dimension of the blob.
Equation (3.2) is the solution for a single blob of practically zero volume fraction. The results diverge greatly from the original values of plasma and blob tensors for blob filling ratios greater than 30 %.
Since the previous descriptions of
$\boldsymbol{e}$
, as in (3.2), are inadequate for our purposes, we seek a more general form of
$e$
that is suitable for a collection of blobs in a plasma. We put forth the following ansatz,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn17.gif?pub-status=live)
It should be noted that the above form of
$\boldsymbol{e}$
, unlike that in (3.2), is not an exact solution. However, the basis for the construction of our ansatz is directly traceable to (3.2). We have assumed that the major axis of the blob is aligned along the externally imposed magnetic field. Subsequently, the electric field is a sum of two terms:
$(1-\unicode[STIX]{x1D70E})$
times the electric field in the background plasma, and
$\unicode[STIX]{x1D70E}$
times the electric field inside the blob.
Care must be given in the definition of
$\unicode[STIX]{x1D64C}_{\text{equiv},B}$
: this can be easily obtained by interchanging definition of plasma and blob in (2.8). This yields:
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn18.gif?pub-status=live)
The indices
$P$
and
$B$
denote, as usual, plasma and blob respectively. Note that the area where
$\unicode[STIX]{x1D64C}_{\text{equiv},B}=0$
is the area where
$\unicode[STIX]{x1D64C}_{\text{equiv},P}\neq 0$
and vice versa.
Our next step is to transform (3.3) into a more usable form. Utilizing (2.5), (2.6) and (2.7) yields:
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn19.gif?pub-status=live)
Where
$H$
denotes the Heaviside step function, valued 1 when the condition in its argument holds, and 0 elsewhere. Fourier transforming the above equation finally gives
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn20.gif?pub-status=live)
where
$\boldsymbol{F}\equiv 4\unicode[STIX]{x03C0}\int ((\sin c|\unicode[STIX]{x1D640}^{-1}\boldsymbol{\cdot }(\boldsymbol{q}-\boldsymbol{q}^{\prime })|-\cos |\unicode[STIX]{x1D640}^{-1}\boldsymbol{\cdot }(\boldsymbol{q}-\boldsymbol{q}^{\prime })|)/|\unicode[STIX]{x1D640}^{-1}\boldsymbol{\cdot }(\boldsymbol{q}-\boldsymbol{q}^{\prime })|^{2})\hat{\boldsymbol{e}}(\boldsymbol{q}^{\prime })\,\text{d}^{3}q^{\prime }$
.
Utilizing (3.1) to express the electric fields as functions of wavenumber and dielectric tensor components, and then integrating (3.6) over the domain of
$\boldsymbol{q}$
provides a system of three equation with 3 unknown variables, which are the components of the dielectric tensor
$\unicode[STIX]{x1D646}_{H}$
. Solving this system provides the necessary solutions for the homogenized tensor.
3.2 Results for the EC (Electron Cyclotron) frequency range
To present our method in contrast to the ones discussed previously, three graphs will take place. This case is for the electron cyclotron plasma frequency
$\unicode[STIX]{x1D714}=170$
GHz. With the components of the plasma and blob matrices given by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn21.gif?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn22.gif?pub-status=live)
Ellipsoidal semiaxes ratios are given by
$a_{y}=0.5$
,
$a_{z}=4$
,
$a_{x}=0.5$
and blob radius
$\unicode[STIX]{x1D70C}=3~\text{cm}$
. The
$z$
axis is along the magnetic field, thus the elongation of the blob toward its direction.
The methods to be compared to ours is the simple mixing formula by Sihvola and an elementary mixing rule given by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_eqn23.gif?pub-status=live)
The results are quite interesting, as it is clear that Sihvola’s method coincides with this simple, non-formal homogenization mixing rule for all components of the tensor. The blue line illustrating our method is interpolated by using values of the filling ratio from 0 to 1. Results from previous methods were not satisfactory in the ‘grey area’ 30–70 % for the original methods of (3.2) with the long wavelength approximation (strong divergence away from the values of the plasma and blob tensors) whereas our method (2.6) behaves smoothly, and oscillates around the other methods, behaving in such a way because of fewer approximations (we incorporate more physics). The oscillation is seen in all six figures 1–6. A dielectric tensor graph of a straight line relation to the filling ratio seems as though it does not incorporate enough physical laws and plasma blob interactions. We propose that our method can be more accurate, and can be improved by using experimental data to verify and interpolate the results into a complete graph.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_fig1g.gif?pub-status=live)
Figure 1. The
$K_{H,\bot }$
component.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_fig2g.gif?pub-status=live)
Figure 2. The
$K_{H,X}$
component.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_fig3g.gif?pub-status=live)
Figure 3. The
$K_{H,\Vert }$
component.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_fig4g.gif?pub-status=live)
Figure 4. Sihvola comparison: the
$K_{H,\bot }$
component.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_fig5g.gif?pub-status=live)
Figure 5. Sihvola comparison: the
$K_{H,X}$
component.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20180905120923689-0046:S0022377818000922:S0022377818000922_fig6g.gif?pub-status=live)
Figure 6. Sihvola comparison: the
$K_{H,\Vert }$
component.
4 Conclusion
We have formulated a modified version of the conventional homogenization theory for describing the permittivity of the turbulent plasma in the scrape-off layer. Using the assumption that the blobs are aligned with respect to the external magnetic field, we proposed the ansatz (3.3) to approximate the total electric field inside the plasma. The formulation is not limited to the long wavelength limit where the wavelength of the RF wave is longer than the scale length of the turbulence. Our model can be extended to higher filling ratios since it has the advantage of incorporating the essential physics of the scrape-off layer and of the RF waves of any frequency. The results show an oscillatory behaviour with varying filling ratio which is quantitatively significant for the off-diagonal terms and sizeable for the diagonal ones. A study on combining the present formalism with the 4 by 4 matrix technique for beam propagation through a scrape-off layer plasma is presently under way. Preliminary aspects of the combined theory have been reported by Valvis et al. (Reference Valvis, Papagiannis, Papadopoulos, Hizanidis, Glytsis, Bairaktaris, Zisis, Tigelis and Ram2017). We intend to apply statistics so as to include filaments and blobs with different shapes, sizes and contrasts, and incorporate data from experiments performed on TORPEX (TORoidal Plasma EXperiment) and TCV (Tokamak à Configuration Variable) (Chellai et al. Reference Chellai, Alberti, Baquero-Ruiz, Furno, Goodman, Manke, Plyushchev, Guidi, Koehn and Maj2018).
Acknowledgements
A.K.R. is supported by the US Department of Energy grant nos. DE-FG02-91ER-54109, DE-FG02-99ER-54525-NSTX and DE-FC02-01ER54648. All the other authors acknowledge partial funding from the EUROfusion Consortium as well as from the National Program for Thermonuclear Fusion of the Association EURATOM/Hellenic Republic.