Measurement and physical interpretation of the mean motion of turbulent density patterns detected by the BES system on MAST

It is necessary to know the true mean velocity of the density patterns to investigate statistical reliability of the cross-correlation time delay (CCTD) method described in section 2.2. Specifically, we investigate how reliable the CCTD method is for different magnitudes of mean velocities and correlation times of the density patterns. We must also evaluate the effect of global modes (i.e., MHD modes) and temporally varying poloidal velocities on the CCTD method. For this purpose, we numerically generate artificial density patterns, random both in space and time, then produce synthetic BES data using the point-spread-functions (PSFs) of the 2D BES system on MAST (as described in Appendix A.2) and compare the inferred flow velocity with the true flow velocity.
We follow a similar approach to the one suggested by Zoletnik et al. [46]. Let the density patterns be described by Gaussian structures both in space and time, namely,

where $x$, $y$ and $t$ denote radial, poloidal and time coordinates, respectively. These numerically generated density patterns are referred to as “eddies” in this paper. Here $N$ is the total number of eddies and the subscript $i$ denotes the $i^{th}$ eddy in the simulation; $\delta n_{0i}$, $x_{0i}$, $y_{0i}$ and $t_{0i}$ are the maximum amplitude and central locations in the $x$, $y$ and $t$ coordinates of the $i^{th}$ eddy, respectively; $\lambda_{x}$, $\lambda_{y}$ and $\tau_{life}$ are the widths of our Gaussian eddies in the $x$, $y$ and $t$ directions; $\tau_{life}$ is the lifetime (or the correlation time) of the eddies in the moving frame; $v_{y}(t)$ is the apparent advection velocity of the eddies in the poloidal direction. Although it is possible to introduce a finite radial velocity shear by making $v_{y}$ a function of $x$, the effect of such shearing rates on the CCTD method is not investigated in this paper, so we will only consider $v_{y}$ that are independent of $x$. The $\cos$ term in the $y$ (poloidal) direction is introduced to model wave-like-structured eddies in the poloidal direction as observed in tokamaks [47]. Note that the envelope (i.e., the $\exp$ term) and the wave structure (i.e., the $\cos$ term) of $\delta n(x,y,t)$ have the same advection velocity $v_{y}(t)$. The central locations of eddies, $x_{0}$, $y_{0}$ and $t_{0}$, are selected from uniformly distributed random numbers, whereas their amplitudes $\delta n_{0}$ are selected from normally distributed random numbers whose standard deviation is one.^§§§It is worth mentioning that there is another scheme of generating such eddies numerically, proposed by Jakubowski et al. [48]. They generated the time series of fluctuating density ($\delta n_{1}$) using the inverse Fourier transform of a broadband Gaussian amplitude distribution in frequency space. Then, a second signal ($\delta n_{2}$) was generated by imposing the desired time-delay fluctuation on the $\delta n_{1}$ such that $\delta n_{2}$ was a time-delayed version of $\delta n_{1}$. This method does not include spatial information for the signals.
The spatial domain of the simulation is $25\>cm$ and $20\>cm$ with the mesh size of $0.5\>cm$ in radial ($x$) and poloidal ($y$) directions, respectively. The time duration of the simulation is $20\>msec$ with a $0.5\>\mu sec$ time step so as to have the same Nyquist frequency as the real 2D BES data from MAST. The widths $\lambda_{x}$ and $\lambda_{y}$ are set so that the full width at half maximum (FWHM) in the radial direction and the wavelength in the poloidal direction are $\sim 8\>cm$ (i.e., $\lambda_{x}=3.53\>cm$) and $\sim 20\>cm$ (i.e., $\lambda_{y}=20.0\>cm$), respectively, which are similar to the measured correlation lengths with the 2D BES system on MAST.^¶¶¶Note that Smith et al. [49] also reported that poloidal correlation lengths of the density patterns are $\sim 20\>cm$ using their 2D BES system on NSTX. The eddy lifetime in the moving frame ($\tau_{life}$) is set to $15\>\mu sec$. However, some of the data sets in this paper have different values of $\tau_{life}$, so the effect of $\tau_{life}$ on the CCTD method can be investigated.
The total number of eddies is $N=20000$. If the eddies are too sparse in the simulation domain, then we may not achieve steady statistical results, while overly dense eddies may cause an effective widening of the specified spatial ($\lambda_{x}$ and $\lambda_{y}$) and temporal ($\tau_{life}$) correlations as many eddies can merge into one larger eddy. Thus, we introduce another control parameter, the spatio-temporal filling factor ($F$), defined as

where $\tau_{ac}$ is the autocorrelation time calculated as

for the generated eddies defined by equation (10). All of our synthetic data was generated so as to $F\sim\mathcal{O}(1)$.
The testing of the CCTD method will involve exploiting what happens if $v_{y}(t)$ has a mean and a temporally varying components. Thus, we generate a temporal structure of $v_{y}$: at each $x$,

where $\left\langle v_{y}\right\rangle$ and $\delta v_{y}$ are the mean and temporally varying velocities, respectively, $\tau_{fluc}$ and $f_{fluc}$ are the lifetime and frequency of $\delta v_{y}(t)$, respectively, and $\tilde{v}_{y}(t)$ is generated from normally distributed random numbers. The RMS (root-mean-square) value of $\delta v_{y}(t)$ denoted as $\delta v_{y}^{RMS}$ will be varied as well as $\left\langle v_{y}\right\rangle$ to investigate the effects of these quantities on the CCTD method. $\tau_{fluc}$ and $f_{fluc}$ allow one to introduce structured temporally varying velocities, while the randomness is kept by $\tilde{v}_{y}$. As one of the causes for the temporal variation of the poloidal velocity is believed to be the existence of geodesic acoustic modes (GAMs)^∥∥∥We do not investigate whether the CCTD method is able to detect such a temporally structured $\delta v_{y}\left(t\right)$ (or GAMs) in this paper, rather we investigate how the existence of these structures affects the CCTD-determined mean velocity. [50], we choose $\tau_{fluc}=500\>\mu sec$ and $f_{fluc}=10\>kHz$ to mimic the GAM features detected by Langmuir probes on MAST [51].
The simulations have been run on a NVIDIA® GeForce GTS 250 GPU card using CUDA programming, which increases the computational speed owing to the highly parallelizable structure of equation (10).

We generate the $i^{th}$ (1 to 8) radial and $j^{th}$ (1 to 4) poloidal channel of the synthetic BES data $I^{ij}\left(t\right)$ by using the calculated point-spread-functions (PSFs) of the actual 2D BES system on MAST [23] and $\delta n\left(x,y,t\right)$ from equation (10) with an additional random noise. Furthermore, a large-scale (in space) coherent (in time) oscillation is included to imitate a global MHD mode. Namely, $I^{ij}\left(t\right)$ is defined as

where $I_{DC}^{ij}$ is the DC value – a typical value of $0.8\>V$ is used for all channels [16]. The rest of the terms are as follows.
$\delta I^{ij}\left(t\right)$ is the fluctuating part of the signal generated from the Gaussian eddies, $\delta n\left(x,y,t\right)$, given by equation (10) and convolved with the PSFs of the 2D BES system :

where $\mathcal{P}^{ij}\left(x,y\right)$ is the PSF of the $i^{th}$ and $j^{th}$ channel of the 2D BES system, normalized so that RMS value of $\delta I^{ij}\left(t\right)$ is $\delta I^{RMS}$. This value is set so that the ratio of $\delta I^{RMS}$ to $I_{DC}^{ij}$ is 0.05. An example of the PSFs for the 32 channels of the 2D BES system on MAST is shown in Figure 9 (white contour lines in the top-left panel).
$I^{ij}_{MHD}\left(t\right)$ models an MHD (global) mode. We assume that the spatial scale of the MHD modes is larger than the BES domain in the poloidal direction, so $I^{ij}_{MHD}\left(t\right)$ does not vary in the poloidal direction. The model MHD signal is generated in a way similar to temporal behaviour of $v_{y}\left(t\right)$ using equation (13), except that the mean value of $I^{ij}_{MHD}\left(t\right)$ is zero and $\tau_{fluc}=250\>\mu sec$. The frequency of the mode $f_{MHD}$ and its RMS value, denoted $I_{MHD}^{RMS}$, will be varied in various tests. The value of $\tau_{fluc}$ here is representative of MHD burst-like fishbone instabilities [52] or chirping modes [53] in tokamaks, for which the spectrum has a finite bandwidth.
$I_{N}^{ij}\left(t\right)$ represents the noise in the signal. As the noise of the 2D BES system on MAST is dominated by the photon noise [19], $I_{N}^{ij}\left(t\right)$ is generated using normally distributed random numbers. Its RMS level is set such that the signal-to-noise ratio ($SNR$) is $300$, which is typical of the 2D BES system on MAST [16].
Figure 8 shows examples of autopower spectra of the synthetic 2D BES data for $\left\langle v_{y}\right\rangle=2.0,5.0,10.0$ and $40.0\>km/s$. The autopower spectrum is calculated as $|FT\{I^{ij}\left(t\right)\}|^{2}$ where $FT\left\{\cdot\right\}$ is the Fourier transform in the time domain.

As defined by equation (2), cross-correlation functions are calculated as time averages of the data. For a $20\>msec$-long synthetic data set containing $N_{total}=40,000$ data points with the sampling time $\Delta t_{sam}=0.5\>\mu sec$, we want to determine $v^{BES}_{y}$ with a time resolution $t_{res}=1\>msec$. First, a cross-correlation function (2) is calculated on a sub-time window of the synthetic 2D BES data containing $N_{\mathcal{C}}$ points, where $N_{\mathcal{C}}<t_{res}/\Delta t_{sam}$. Then, such cross-correlation functions are averaged over $N_{avg}$ consecutive sub-time windows where $N_{avg}=(t_{res}/\Delta t_{sam})/N_{\mathcal{C}}$ so that an averaged cross-correlation function is obtained at every $t_{res}$. In this paper, we use $N_{\mathcal{C}}=80$, so $N_{avg}=25$.
Denoting $f(t)$ and $g(t)$ the time series over a sub-time window from two poloidally separated synthetic 2D BES channels, the cross-correlation function (2) for this sub-time window is:

for any integer $r$ with $|r|<N_{\mathcal{C}}-1$. Finally, by averaging $\mathcal{C}_{sub}$ for $N_{avg}$ consecutive sub-time windows we obtain the smoothed averaged cross-correlation function $\mathcal{C}(r\Delta t_{sam})$ from $1\>msec$-long data points.
The CCTD method has a serious limitation due to the fact that the sampling time $\Delta t_{sam}$ is finite. In order to calculate $v_{y}^{BES}$ using only two poloidally separated channels, a line is fitted through two points on a $\left(\Delta y,\tau_{peak}^{cc}\right)$ plane as shown in Figure 1(b). The first point is located at $\left(\Delta y,\tau_{peak}^{cc}\right)=\left(0,0\right)$ by definition, and the second point at $\left(\Delta y,r\>\Delta t_{sam}\right)$. Then, possible values of $v_{y}^{BES}$ are restricted to $\Delta y/\left(r\Delta t_{sam}\right)$ where $r$ is an integer. For the 2D BES system on MAST, using two adjacent poloidal channels ($\Delta y=2.0\>cm$) with a sampling time $\Delta t_{sam}=0.5\>\mu sec$, the possible values of $v_{y}^{BES}$ are limited to $40.0,\>20.0,\>13.3,\ldots\>km/s$ for $r=1,\>2,\>3,\ldots$. Such a limitation may be mitigated by using four poloidally separated channels. However, using four channels is not always possible if the channels that are farthest apart are not correlated. To resolve this issue, we use a second-order polynomial fit on the cross-correlation function $\mathcal{C}(r\Delta t_{sam})$ to locate its global maximum: if $r_{peak}$ is the point where the discrete cross-correlation function $\mathcal{C}(r\Delta t_{sam})$ is maximum, we use the three values of $\mathcal{C}(r\Delta t_{sam})$ at $r=r_{peak}$, $r_{peak}-1$ and $r_{peak}+1$ to fit a second-order polynomial. The “true” maximum is found from this fit. We denote the time delay at which this maximum is reached by $\tau_{peak}^{cc}$.

For a given set of $20\>msec$-long synthetic 2D BES data, we calculate $v_{y}^{BES}$ with the time resolution of $1\>msec$ (Appendix B.1). Furthermore, we do this at three different radial locations^******As described in Appendix A, $v_{y}\left(t\right)$ are identical at all radial locations. One column in the middle and two columns from the edges of the 2D BES channels are used. so that the average of $v_{y}^{BES}$, denoted $\left\langle v_{y}^{BES}\right\rangle$, can be calculated using $60$ values of $v_{y}^{BES}$. To make quantitative comparisons between $\left\langle v_{y}^{BES}\right\rangle$ and $\left\langle v_{y}\right\rangle$ defined in equation (13), we define four types of error.
The normalized bias error

is a quantitative measurement of the systematic discrepancy between the measured and the true value. The normalized random error

quantifies the degree of fluctuation in the measured $v_{y}^{BES}$ with respect to $\left\langle v_{y}^{BES}\right\rangle$. This value may depend on the MHD contribution in equation (14) and the temporally varying poloidal velocity $\delta v_{y}(t)$ in equation (13).
Furthermore, as linear fitting is done to determine $v_{y}^{BES}$ (see Figure 1), two other types of error are present. The slope of a linear fit can be denoted as $v_{y}^{BES}\pm\delta v_{fit}$ where $\delta v_{fit}$ is a degree of the uncertainty of the least-square fit^††††††In Figure 1, we plotted $\tau_{peak}^{cc}$ as a function of $\Delta y$ and determined $v_{y}^{BES}$ as the inverse of the slope of a fitted line. Operationally, we actually plot $\Delta y$ as a function of $\tau_{peak}^{cc}$ so the slope of a fitted line is the $v_{y}^{BES}$.. Then, the normalized mean of $\delta v_{fit}$ is

and the normalized random error in $\delta v_{fit}$ is

These two uncertainties together provide an estimation of how well a linear line is fitted to given data points. For example, if the assumption that $\tau_{life}$ is long enough so that all four poloidally separated channels observe the same eddies is not satisfied, then $\hat{\sigma}_{mean}^{fit}$ becomes large. On the other hand, if this assumption is occasionally satisfied, then $\hat{\sigma}_{rand}^{fit}$ exhibits such events because $\delta v_{fit}$ will then be small compared to its average. Note that error bars of the CCTD-determined apparent velocities in Figures 4 - 7 show $\left\langle\delta v_{fit}\right\rangle$.
In the following sections, these four types of error will be evaluated for various values of $\left\langle v_{y}\right\rangle$ and $\tau_{life}$, and various ranges of $I_{MHD}^{RMS}$, $f_{MHD}$ and $\delta v_{y}^{RMS}$.

To investigate the reliability of the CCTD method described in Appendix B.1 for estimating $v_{y}^{BES}$, we generate a number of synthetic 2D BES data sets with various values of $\left\langle v_{y}\right\rangle$ while keeping all the other parameters in equations (10), (13), (14) and (15) constant. In real experiments, there is almost always some temporal variation of $v_{y}$, thus the RMS value of $\delta v_{y}$ in equation (13) is set to $5\>\%$ of $\left\langle v_{y}\right\rangle$ in this subsection. The synthetic 2D BES data are frequency-filtered to suppress the noise before the cross-correlation functions are calculated. Figure 10 shows examples of (a) $v_{y}\left(t\right)$ generated according to equation (13) with $\left\langle v_{y}\right\rangle=10.0\>km/s$ and (b) the original (black) and frequency-filtered (red) autopower spectra of a generated synthetic signal. Here, the noise cut-off level is set to be the $5$ times the averaged autopower level above $900\>kHz$ (green dashed line).

(1) For $\left\langle v_{y}\right\rangle\lesssim 5.0\>km/s$, the CCTD method is not reliable. This is due to the fact that eddies do not live long enough to be detected by all the poloidally separated channels. Indeed, it was a priori clear that $\left\langle v_{y}\right\rangle<\Delta y/\tau_{life}$ could not be measured. This translates to $\left\langle v_{y}\right\rangle<4.0\>km/s$ for $\Delta y=6.0\>cm$ and $\tau_{life}=15.0\>\mu sec$, so our results are consistent with this simple criterion.
(2) The CCTD method usually overestimates $\left\langle v_{y}\right\rangle$ (i.e., $\hat{\sigma}_{bias}>0$). This can be explained by the effective channel separation distance ($\Delta y$) being in fact slightly less than $2.0\>cm$ because of the overlapping of the PSFs, as shown in Figure 9.
(3) The limitation of the CCTD method due to the finite $\Delta t_{sam}$ is successfully overcome by fitting a second order polynomial to the cross-correlation function, as explained in Appendix B.1.

As explained in section 2.2.2, the CCTD method for determining $\left\langle v_{y}\right\rangle$ is based on the idea that the peak of the cross-correlation function occurs at $\tau_{peak}^{cc}=\tau_{prop}$, where $\tau_{prop}=\Delta y/\left\langle v_{y}\right\rangle$ is the propagation time of the fluctuating density patterns between detectors poloidally separated by the distance $\Delta y$. However, $\tau_{peak}^{cc}$ will not coincide with $\tau_{prop}$ if the lifetime $\tau_{life}$ of the fluctuations is not long compared to $\tau_{prop}$. The failure of the method for $\left\langle v_{y}\right\rangle<5.0\>km/s$ illustrated in Figure 11 is an example of what happens when $\tau_{prop}$ is too large. Here, we investigate the effect of $\tau_{life}$ on $\tau_{peak}^{cc}$ quantitatively, via a systematic $\tau_{life}$ scan of the synthetic BES data.
Two values $\left\langle v_{y}\right\rangle=5.0$ and $20.0\>km/s$ are chosen for this study. For $\left\langle v_{y}\right\rangle=5.0\>km/s$, $\tau_{prop}=4.0$, $8.0$ and $12.0\>\mu sec$ with $\Delta y=2.0$, $4.0$ and $6.0\>cm$, respectively; for $\left\langle v_{y}\right\rangle=20.0\>km/s$, they are $1.0$, $2.0$ and $3.0\>\mu sec$. The peak time $\tau_{peak}^{cc}$ is found using the polynomial fitting method described in Appendix B.1, and $\left(\tau_{prop}-\tau_{peak}^{cc}\right)/\tau_{prop}$ as a function of $\tau_{life}$ is plotted for three different values of $\Delta y$ in Figure 12. It shows that $\tau_{peak}^{cc}$ underestimates the true $\tau_{prop}$ for small values of $\tau_{life}$, leading to an overestimation of the $\left\langle v_{y}\right\rangle$, consistent with the results shown in Figure 11. It is encouraging, however, that even relatively low velocities of just a few $km/s$ can be determined by the CCTD method with reasonable accuracy ($\sim\>20\%$).