Magnetized plasma pressure filaments: Analysis of chaotic and intermittent transport events driven by drift-Alfvén modes

Filamentary structures are a persistent feature in many magnetized plasma environments such as laboratory, solar atmospheric, planetary magnetospheres, and astrophysical plasmas. Recent experiments in magnetized dusty plasmas have observed filamentary structures between electrodes that are aligned with an external magnetic field and disrupt the uniformity of the plasma.¹ Radially propagating filamentary structures have been observed in the boundary layer regions of fusion plasmas and are believed to significantly enhance the transport across the scrape-off layer (SOL).² In the last two decades, there has been mounting experimental evidence from toroidal and linear magnetized plasma devices^3–5 that coherent structures in the plasma turbulence, referred to as blobs or magnetic field-aligned filaments, lead to intermittent convective cross field transport of particles, energy, and momentum.⁶ 

Intermittent fluctuations that are frequently observed in laboratory and space plasmas are typically described as sporadic large amplitude bursts within a quiescent background that lead to a non-Gaussian probability distribution function (PDF) in turbulence measurements. These fluctuation measurements are made using plasma diagnostics, such as Langmuir probes,⁷ beam-emission spectroscopy (BES),⁸ and gas puff imaging.⁹ Intermittency is often associated with turbulent coherent structures, an example of which is the drift-Alfvén vortex filament in a non-uniform and strongly magnetized plasma.⁶ A deterministic model based on chaotic advection due to the interaction of drift-Alfvén waves can generate a burst train of Lorentzian pulses that leads to an exponential frequency power spectral density.¹⁰ An alternative approach uses a stochastic model consisting of self-similar pulses with exponential profiles arriving at random times.¹¹ Using the advanced tools of time series analysis, such as the permutation entropy and the entropy–complexity (CH-plane), it has been possible to distinguish between these two dynamical systems.¹² 

In basic studies of electron heat transport, electron temperature filaments have previously been investigated under controlled conditions using a single isolated heat source embedded in a large linear plasma device.^13–16 In these experiments, a low-voltage electron beam was injected into a strongly magnetized, colder, afterglow plasma. This produced a long (∼8 m) and narrow (∼10 mm diameter) filament of elevated temperature (∼20 times the background), isolated from the walls of the chamber (1 m diameter). The experiments established that there is a transition from a period of classical transport¹⁷ (due to Coulomb Collisions) to one of anomalous transport.¹⁸ In this latter phase, localized drift-Alfvén eigenmodes were driven unstable by the temperature gradient in the filament edge. As the plasma conditions changed, the highly coherent eigenmodes evolved into broadband drift-Alfvénic turbulence.^10,19 These studies were extended to include several interacting electron temperature filaments that exhibited stronger nonlinear behavior.²⁰ 

The purpose of this work is to further investigate the properties of nonlinear convective transport due to drift-Alfvén waves in controlled multi-filament plasma structures. We perform the statistical analysis of the fluctuations and generate the PDF, thus revealing the emergence of intermittency in the system. The intermittent events are determined to have a temporal Lorentzian shape and the many thousands of such events populating the time series are characterized according to width, location, and amplitude. By making use of cross-conditional averaging of pulse events using a reference probe, we have characterized the spatiotemporal evolution of the events as localized structures that propagate radially and azimuthally away from the filaments. The pulse analysis and spatial reconstruction reveal that interaction among drift-Alfvén modes is correlated with the event occurrence. Last, we use the complexity–entropy plane (CH-plane) analysis to discriminate between stochastic and chaotic signals and show that the intermittent events are driven by deterministic chaotic dynamics.

The paper is organized as follows: Sec. II describes the experimental setup and measurement methods. Section III presents the experiment results, including a description of the global drift-Alfvén eigenmodes that are the source of the turbulence, a statistical analysis of the intermittent pulse events, and a study of deterministic chaos in the filament system; and Sec. IV is a summary and discussion of the results presented.

The experiment is performed in the upgraded Large Plasma Device (LAPD)²¹ at the Basic Plasma Science Facility (BaPSF) at the University of California, Los Angeles (UCLA). The LAPD, shown schematically in Fig. 1(a), is a linear plasma device that produces a magnetically confined cylindrical plasma column that is 18 m long and 60 cm in diameter. The cylindrical chamber is encircled by electromagnets that produce the background magnetic field, fixed at $B 0 = 0.1$ T (1000 G) for the experiments herein. The He gas filling the chamber has a fill pressure of $13.3 × 10 − 5$ Torr (⁠ $1.77 × 10 − 5$ kPa). At one end of the chamber is the main plasma source in the form of a hot barium oxide (BaO) cathode and a mesh anode 50 cm away axially; a bias of −65 V is applied at the cathode, with respect to the anode, and thermionic electrons accelerate into the main plasma chamber and collisionally ionize the He gas to form the plasma column. The main plasma discharge lasts for 12 ms and has a typical electron temperature of $T e ∼ 5$ eV, density of $n ∼ 2 × 10 12$ cm^{− 3}, and ion temperature, T_i, less than 1 eV. Once the main discharge bias is eliminated, the plasma transitions to the afterglow phase where the electron temperature rapidly cools to T_e around 0.25 eV, whereas the density decays exponentially throughout the afterglow with a time constant on the order of 10 ms. The LAPD is a pulsed device and the process of making the plasma is repeated at a 1 Hz repetition rate; the plasma conditions are highly reproducible shot-to-shot.

At a distance of 1500 cm away from the BaO cathode, three 3-mm-diameter crystal cathodes of cerium hexaboride (CeB₆) mounted on probe shafts are inserted into the plasma. Figure 1(d) shows one of the mounted crystals next to an American dime for scale. The crystals are supported by current carrying wires mounted on a ceramic base that has been slightly modified to a semi-circle shape. The wire mounts are connected to insulated wires inside the probe shaft and the outputs are accessible from the end of the probe outside the LAPD vacuum chamber. CeB₆ has a low work function and high electron emissivity when it is heated to an operating temperature of around 1400 °C; the heating is accomplished ohmically with around 10 W of power to the crystal. To align all of the crystals in the LAPD, one probe is inserted from each side of the device and one from the top shown schematically in Fig. 1(a); the axial location of the filaments is taken as $z 0 = 0$ cm. The crystals can be arranged arbitrary close together in the (x, y) plane; in addition to the modifications to the crystal bases, this was facilitated by staggering the probes by a few cm along the z direction and slightly rotating them as seen in Fig. 1(b). The image in Fig. 1(c) taken through a surveying scope aligned along the axis of the LAPD shows the crystals aligned in a quasi-symmetric filament bundle with diagnostic probes visible in the distance. Each crystal has independent heating and biasing circuitry.

A negative bias is applied between each crystal and the mesh anode starting 2 ms into the afterglow phase; the time of application of the bias is taken as t = 0 and the bias last for 20 ms. When the bias is applied between the hot crystal cathodes and the mesh anode at the other end of the device, electrons are emitted from the CeB₆ crystal. The potential difference between the crystals and anode is kept below 20 V to ensure the energy of the thermionic electrons is below the ionization energy of helium. The emitted electrons thermalize in the afterglow plasma after a few mean free paths creating a heated region less than 1 m in extent and a few mm in diameter. The asymmetry in the perpendicular and parallel thermal transport coefficients in a magnetized plasma causes the heated regions to rapidly form filamentary structures of elevated temperature less than 10 m in length with symmetric Gaussian-like transverse profiles about 1–2 cm in diameter. The peak temperature near the heat sources is 3–5 eV depending on the bias voltage and decreases continuously toward the end and edges of the filaments where the temperature equilibrates with the cold background plasma.

Access ports are spaced 32 cm apart axially along the length of the LAPD chamber. The measurements presented in this paper are all collected using small Langmuir probes inserted through the ports and biased to sample the temporal evolution of the current drawn from the plasma. Since the LAPD has a high repetition rate and high reproducibility a probe located at an axial position z can be placed at a position in the transverse (x, y) plane, collect several nearly identical shots, and then be moved to a new position in the plane to repeat the process; this allows the collection of large 2D data planes from an ensemble of plasma shots. The process of moving the probes and collecting the data is entirely automated with the probes mounted on probe drives capable of less than 1 mm accuracy. In these experiments, measurement probes for collecting 2D planes of data are inserted at distances $z 1 = 256$ cm and $z 2 = 544$ from the crystal cathodes, as seen in 1(a). Additionally, two reference probes are inserted manually at z_R = 320 cm, one from each side of the device, and left in the same positions throughout the experiment to use as a reference for correlating the plasma shots in the moving planes.

The Langmuir probes are used to collect two types of measurements, ion saturation current, and temporal Langmuir sweeps. Ion saturation current, $I sat$⁠, is collected at the probe face when the probe is biased well below the plasma potential; the measured current is proportional to the plasma density and square root of the electron temperature, i.e., $I sat ∝ n T e$⁠.²² The $I sat$ measurements are often decomposed into fluctuating and time-averaged components,

using digital filtering to remove frequency content below 1 kHz from the $I sat$ signal to get $δ I sat$⁠. The brackets indicating time-averaging, $⟨ ⋯ ⟩$⁠, are dropped for convenience where fractional amplitudes are used, that is $δ I sat / I sat ≡ δ I sat / ⟨ I sat ⟩$⁠. Temporal Langmuir sweeps are used to obtain measurements of the electron temperature, plasma density, and space potential, Vs, from a characteristic I–V curve.²³ The probe voltage is swept across a voltage range in a continuous sawtooth pattern with a period typically in the range of 200–400 μs. This method is useful for acquiring information about the time-averaged plasma parameters but has poor temporal resolution and cannot yield information about fast fluctuations in the plasma.

In this section, the results from the experiments are presented, beginning with a characterization of the dominant fluctuations driven by the cross field pressure gradients in the magnetized thermal filaments. Using Langmuir probes, transverse planes of $I sat$ signals taken at various axial locations are used to map out the global eigenmode structure and relate these fluctuations to the gradients. This is a prelude to the statistical analysis of the fluctuations beginning with an analysis of the amplitude probability distribution (PDF), which can assess the degree of intermittency in the turbulence.

The ability to independently control the location and strength of heat sources from the crystal cathodes described in Sec. II has allowed us to investigate the interaction in close proximity. Starting with a single source, we have studied the interaction with a second source of equal strength and then with a third source placed in such a way that a triangular pattern was created. The $I sat$ planes located at axial position z₁, for each of these cases is shown in Fig. 2. In this figure, a dual spatial scale is added to indicate the separation in physical units (cm) and normalized units (electron skin depth, δ_e). It is interesting to note that for filaments of ∼1 cm width (∼2δ_e) there is relatively strong interaction at separations of ∼3 cm (⁠ $∼ 5 δ e$⁠), or greater than twice the filament widths, in both the two and three filament arrangements. This interaction is characterized by the formation of tail-like structures extending from each filament in Figs. 2(b) and 2(c) driven by $E × B$ flows and asymmetric drift wave modes that then develop on the tails.²⁰ For the remainder of this paper, we will focus on the tri-filament configuration where the distance of the cathodes from a central origin is adjusted to be ∼5 mm or in terms of filament separation, approximately $1.5 δ e$⁠.

As a first step, we compare the global profiles of electron temperature, density, and space potential for the case of a single filament and the three filaments in close proximity. In the first case, we take an azimuthal average at different radial positions from the center of the single filament and in the latter case, we take an azimuthal average from the center of the triangular pattern. It has previously been demonstrated that the azimuthal average of the filament bundle is a good approximation due to the quasi-symmetry of the tight configuration.²⁰ A comparison of the profiles for the two cases is shown in Fig. 3. The electron temperature and space potential of the multi-filament configuration in Fig. 3(b) are broadened as expected with the temperature peaking off-axis and exhibiting a gradient reversal within the core region. In both cases, the observed well in the space potential of the plasma and associated radial electric field give rise to sheared $E × B$ flows in the azimuthal direction; the potential profiles are a result of the CeB₆ crystals and a recent model for the emissive cathodes²⁴ accurately describes the observed profiles.²⁰ An interesting feature is the reversal of the density gradient in the tri-filament configuration, which peaks outside the higher temperature region of the bundle as compared with the single filament density that peaks near the center. In both cases, the density is relatively flat, thereby indicating it is the electron temperature gradient that is the dominant free energy source for the drift-Alfvén mode instability which is the main source of the fluctuations. For our further analysis, it also implies that $I sat$ and its changes are dominated by the electron temperature.

In the three filament case, Fig. 4(a) presents the time traces of the ion saturation current, $I sat$⁠, shown at several locations within the filaments and these are marked in the $I sat$ planes taken at two different times and two axial locations in the lower panels. It is clear from the $I sat$ color maps Figs. 4(b) and 4(d) that at early time, t = 3 ms, the transverse profile of the heating is non-uniform axially with an apparent merging of the filaments at the further port, z = z₂, from the source. Starting at around t = 12 ms the merged bundle at z₂ begins to broaden, starts to separate at t = 13 ms, and by t = 14 ms three distinct filaments are evident at both ports, shown in Figs. 4(c) and 4(e).

Using the time window of $t = 3 − 6$ ms, the power spectra of the fluctuations $δ I sat$ are displayed in Figs. 5(a) and 5(b) for the two different axial planes at $z 1 = 256$ and $z 2 = 544$ cm, respectively. The radial locations for the single position measurements and azimuthal averages are in the vicinity of the maximum thermal gradient at each axial location. Several dominant frequency peaks are evident along with a background exponential frequency spectrum. The origin of the broad exponential spectrum is discussed in Sec. III C. The spectral peaks, indicated by the dotted black lines, are chosen as the frequencies at which to compute cross-correlations between the moving probe at a particular plane and a nearby reference probe at z_R, thus illuminating the mode structures at each frequency in the transverse plane.

The results of this analysis are shown in Figs. 6(a), 6(f), and 6(k) for the axial location z₁ and Figs. 7(a), 7(f), and 7(k) for z₂. It was previously shown that these modes correspond to drift-Alfvén waves driven by the outer gradient of the filament bundle; a 1D symmetric stability analysis using the temperature, density, and shear flow profiles was used to accurately reproduce the mode structures and predict the frequencies.²⁰ Here, we perform a Fourier decomposition of the quasi-symmetric mode patterns into azimuthal mode numbers at each frequency; the four most dominant modes are shown in the panels to the right of each cross correlation. In Fig. 6, the decomposition of the three dominant spectral peaks at 3, 6.5, and 11 kHz reveals the dominance of mode numbers m = 1, m = 2, and m = 3, respectively. However, the observed modes are in fact superpositions of several coherent modes; for example, the m = 2 mode is nearly as prevalent as the m = 1 at 3 kHz and the pattern at 14 kHz is a superposition of comparable modes at m = 2, 3, and 4. This is evidence of some form of nonlinear interaction between the drift-Alfvén modes and is a common route to turbulence in magnetized plasmas.^25,26 Figure 7 shows that from the source at z₂ the dominant power resides in the m = 1, m = 3, and m = 2 eigenmodes for each of the dominant frequencies 6.5, 14, and 28 kHz, respectively. This shift in azimuthal mode dominance is correlated with axial inhomogeneity. Again, while there is a dominant mode at each frequency, there is comparable power in other mode numbers that together form the observed mode patterns.

In Fig. 8, the radial structure of the largest decomposed modes at each frequency are shown along with the azimuthally averaged and normalized $I sat$ profiles at z₁ and z₂. As the off-axis peak in the profile at z₁ shifts to the center at z₂, the peaks in the mode amplitude also move inwards and track the steepest region of the outer gradient of the bundle. The nonlinear mode coupling observed at z₁ and z₂ provides a mechanism for the turbulent fluctuations in $I sat$⁠, observed in Fig. 4(a) and discussed in Sec. III B, that are characterized by intermittent events localized to the outer gradient of the filamentary structures. Here, we note that the inner gradient (r < 0.4 cm) at z₁ observed in Fig. 8(a) also shows some evidence of a low amplitude, higher frequency wave mode in the range of 25–30 kHz. The cross correlation at these frequencies (between z_R and z₁) is incoherent and suggests the mode is localized to the center of the bundle and the reference probes on the outer gradient do not pickup the signal. It can be surmised that this is a low-amplitude drift-Alfvén wave driven by the highly asymmetric inner gradient; however, we leave discussion and analysis of this mode to future experiments that could place a reference probe into the inner bundle to confirm its existence.

Intermittency can be identified by non-Gaussian amplitude probability distribution functions (PDFs) of the measured signals. An increase in the frequency of large intermittent events within a signal will significantly impact the tails of the distribution. A larger negative tail in the distribution indicates depletion (or “hole”³) events and a larger positive tail indicates enhancement events; the depletions or enhancements in $I sat$ could be density or energy. The color plot in Fig. 9(a) shows the azimuthally averaged amplitude PDF on a logarithmic scale at each radius of the filamentary bundle at z₁ between t = 3 ms and t = 6 ms with 10 shots at each position in the (x, y) plane. The fluctuation amplitudes $δ I sat / I sat$ are normalized to the root mean square deviation of the distribution, σ. Shown in the left panels, in Figs. 9(b)–9(e), are temporal traces of the fluctuations in single shots at radii r = 0.1, 0.5, 1.2, and 2.0 cm and indicated by dotted lines in Fig. 9(a). The right panels, in Figs. 9(f)–9(i), display the PDF of the normalized amplitude for the same radii; the PDF in Fig. 9(i) corresponds to a cumulative PDF of the hatched region (⁠ $r ≥ 2.0$⁠) in Fig. 9(a) where the fluctuation amplitude PDFs are all similar. It is important to note that while the filament bundle at z₁ has some azimuthal asymmetry in the $I sat$ profile, the azimuthal variation of the PDF was found to be minimal.

Near the center of the multi-filament structure, Figs. 9(b) and 9(f), the signal is characterized by nearly Gaussian fluctuation amplitudes with a small enhancement of the negative tail, indicating some tendency for depletion events. The 2D PDF shows that the frequency and size of these depletion events grows further from the center and peaks at around r = 0.5 cm; this is also where the radial $I sat$ profile peaks in Fig. 8(a). The PDF at r = 0.5 cm is more clearly seen in Fig. 9(g); in the corresponding time trace in Fig. 9(c), it is observed the events are short-lived in time only lasting on the order of several μs.

Beyond $r ≈ 0.8$ cm, the PDF switches to an asymmetry skewed toward positive events on the outer gradient. At r = 2.0 cm; and beyond the PDF is well described as Gaussian as evidenced by the fluctuations at 2.0 cm in Fig. 9(e) and the spatial cumulative PDF in Fig. 9(i). The extreme of the positive skewness to the PDF occurs at r = 1.2 cm and corresponds to the amplitude peak of the drift-Alfvén modes in Fig. 8(a). The single shot at r = 1.2 cm in Fig. 9(d) and associated PDF in Fig. 9(h) shows the positive events are large pulses with widths on the order of tens of μs. A similar PDF analysis is carried out at z₂ and is shown in Fig. 10; the format is the same as Fig. 9. Here, the fluctuations are nearly Gaussian at the center of the combined filament structure. A sharp transition to positive intermittent transport events occurs on the gradient between r = 0.4 cm and r = 1.0 cm, with the extreme of the skewness again coinciding with the peak in the drift-Alfvén amplitudes. The fluctuations become low amplitude and nearly Gaussian beyond r = 1.5 cm, similar to beyond r = 2.0 cm at z₁.

At both axial locations, z₁ and z₂, large positive intermittent events are located on the outer gradient of the filament bundle where there are drift-Alfvén modes and the existence of a broadband exponential frequency spectra, as seen in Fig. 5. Previous work has highlighted the observation of Lorentzian-shaped pulses in the turbulent regime created by only a single filament.^10,19 These pulses were observed on the gradient of the single filament and are linked to the observation of anomalous cross field energy transport rates, exponential frequency spectra, and deterministic chaos.^12,27–29 Thus, a closer inspection of the positive events on the outer gradients at z₁ and z₂ is warranted.

In this section, we present a statistical analysis of these transport events and use conditional averaging to perform a spatial reconstruction. Identification of the positive intermittent pulses embedded in the present data is straightforward due to the large amplitudes above the statistical variation. A search algorithm is used to identify pulses in time series of length 8192 points with a time step of 1.28 μs across the time range t = 3 to $∼ 13$ ms; the time is restricted to the range where the filaments remain merged at z₂. The time signals are first low pass filtered to remove frequency content above 35 kHz. A dynamic threshold amplitude is used to isolate maxima corresponding to the large intermittent events; a threshold of twice the mean of the absolute signal amplitude over a moving 1 ms window was found to be suitable. There is no doubt that some pulses are not identified—some may be small and embedded within the statistical variation and others obscured or distorted by wave fluctuations in their temporal vicinity—this is an unavoidable outcome and many thousands of pulses are still identified within all of the shots.

Lorentzian-shaped pulses can be distinguished from Gaussian or other similarly shaped events by the presence of an exponential frequency spectrum associated with the pulses. A general form for pulsed events is described by the stable distribution given by the following function in the Fourier domain,³⁰ 

where

In Eq. (2), t₀ is a location parameter, τ is a scale parameter (⁠ $τ > 0$⁠), β is a skewness parameter (⁠ $− 1 ≤ β ≤ 1$⁠), and α is known as the shape or stability parameter (⁠ $0 < α ≤ 2$⁠). For zero skewness (i.e., β = 0), the distribution is symmetric. The distribution is equivalent to a Levi distribution if $α = 0.5$⁠, a Gaussian distribution if α = 2, and a Lorentz (or Cauchy) distribution if α = 1. The shape parameter, α, is important for the discussion of the form of the function in the Fourier domain; only when α = 1 is the magnitude of the function, $| L |$⁠, exponential in frequency and this remains so for any value of β. The exponential decay rate, or the slope on a log-linear spectrum, is proportional to the parameter τ. Thus, an exponential frequency spectrum can be evidence of Lorentzian-shaped signals in the time domain, and vice versa. Maggs and Morales²⁷ have argued that the log –log formats of power spectra have obscured the exponential frequency spectrum of intermittent turbulence observed in several magnetic confinement plasma devices; additionally, the intermittent pulses in these devices, including similar filamentary experiments in the LAPD, are Lorentzian shaped²⁸ and the source of the exponential frequency spectra.¹² 

The broadband exponential spectra in Fig. 5 are strong evidence that the pulses identified herein are also of Lorentzian form. In Fig. 11(a), an example of a probe signal with a single intermittent event at z₁ and r = 1.5 cm is displayed along with the filtered time trace used for identifying the pulses systematically. Figure 11(b) shows the same filtered pulse with a fit to Eq. (2) using α = 1 and additional amplitude and offset parameters. The pulse is accurately described by a Lorentzian function with $τ = 15.3 ± 0.4$ μs and a positive skewness of $β = 0.53 ± 0.03$⁠; the error is described by 95% confidence intervals. The choice of fixing α = 1 is justified by the presence of the exponential spectra, however, we note that performing the fit with a variable α results in similar τ and β with $α = 1.01 ± 0.04$ and further confirms the Lorentzian shape of the pulses. Though not shown here, the resulting fit in Fig. 11(b) matches the exponential spectrum in the frequency domain and a comparison will be made shortly.

Verifying an accurate fitting of each pulse to Eq. (2) is not practical and each pulse can be quickly characterized by a normalized amplitude (⁠ $δ I sat / I sat$⁠) and a full width at half maximum (FWHM), $Δ t FWHM$⁠. The Lorentzian scale parameter τ is related to $Δ t FWHM$ by the approximation $Δ t FWHM ≈ 2 τ$⁠, where the approximation becomes an equivalence for a symmetric pulse with β = 0. Figures 12 and 13 show the pulse amplitude PDFs and bivariate PDFs of pulse width and amplitude at different axial locations and radii.

In Fig. 12(a), we can see that the normalized amplitude PDFs are well described by lognormal distributions with decreasing mean amplitudes further from the outer gradient. For $r ≈ 1.8$ cm the pulses identified on the reference probe at z_R are used instead of those on the moving probe at z₁ since the larger number of pulses identified provides a better indication that the lognormal fit is appropriate; the slight discrepancy for $δ I sat / I sat ≲ 0.1$ is likely due to the pulses being too small to easily identify. In the bivariate distribution of Fig. 12(b), the pulse width is narrowly distributed at r = 1.0 cm (note the logarithmic PDF scale) and the distribution spreads out further down the gradient in Figs. 12(c) and 12(d). At z₂ in Fig. 13(a) similar lognormal distributions in the amplitude are observed; the amplitude distributions again shift to lower mean amplitudes further from the center and become more narrow on the edge of the filamentary structure. The bivariate distributions in Figs. 12(b)–12(d) show less evidence of the pulse width distribution spreading out further from the center; the filaments remain merged at this axial location and the tighter configuration may be the cause of the spatial homogeneity of the event widths.

The pulses identified can be used to construct the average pulse shape; Figs. 14 and 15 show this result for a radial location on the gradient at z₁ and z₂. In both figures, the top panels (a) show the time trace of the average pulse fit to Eq. (2) and the bottom panels (b) show the pulse and fit in the Fourier domain. The inset figures in the top panels show a larger temporal extent where outside the main pulse the averaging of the shots results in a nearly flat trace. However, on average, in both cases the pulses are within a noticeable depression in $I sat$⁠. In the power spectra, the exponential character is clear and the fits are an obvious match to the slope. Also present in the spectra is evidence of the peaks at 6 and 14 kHz, corresponding to the drift-Alfvén modes at both axial locations, and suggests the wave modes may be correlated with the pulse generation.

The placement of the reference probes on the gradient allows for conditional averaging of the temporal signals on the moving probe at z₁. When a pulse is detected on one of the reference probes, the location of the peak is used as a temporal marker to which the corresponding signal on the moving probe is referenced to as t = 0. This is repeated for each plasma shot at all locations of the moving probe. In this way, a spatial reconstruction of the average pulse event at z₁ is made; the result is shown in Fig. 16. It is important to note that in this analysis, the observed spatial structure of the pulse is a reconstruction from many pulse events detected at the reference probe that do not necessarily represent identical spatially situated events, either in the transverse (x, y) plane or axially in z. Some of the observed events on the reference probe could be the end or beginning of a pulse, near the edge of a pulse in the transverse plane or axially, and some pulses may not even extend to z₁ and may wash out or distort the reconstruction. To help mitigate this issue, only the largest 50% of pulses detected at z_R were used as these pulses are more likely to represent events with maxima localized to the vicinity of the reference probe; using all of the pulses does result in a similar reconstruction, albeit one that is diminished in amplitude and more spread out in the transverse plane, as expected. Additionally, given that a spatial reconstruction at z₁ is possible when the events occurred at z_R, it can be concluded that the axial extent of the events is at least $Δ z 1 , R = z R − z 1 = 64$ cm. Performing the same conditional averaging at z₂ results in no coherent reconstruction of the pulse event, indicating the pulses at z_R are not correlated with pulses at z₂, and an upper bound on the axial extent of the blob structures can be set at $Δ z R , 2 = z 2 − z R = 224$ cm.

In Fig. 16(a), the time traces track the path of the reconstructed event and each has been fit to Eq. (2) showing the persistence of the Lorentzian pulse shape throughout the evolution. The panels of Figs. 16(b)–16(e) show frames of the normalized fluctuation level with the event indicated by the dashed contour line and panels (f)–(i) show the normalized ion saturation current with the same contouring around the event structure; the colored crosses in each frame indicate the location of the time traces in (a) and the dashed lines in (a) show the time of each frame. A clear feature of the fluctuation panels is that the phase of the drift-Alfvén waves is also correlated with pulses occurring at z_R. At time t₁, there is an obvious mixed mode structure between a lower mode number of m = 2 and a higher mode number of perhaps m = 3 or 4. The pulse event itself evolves out of an arm of the drift wave and grows to a large amplitude before dissipating. At its largest extent in the transverse plane [Fig. 16(d)], the event is on the order of $2 δ e$ across. In fluctuation panels for t₂ and t₃ the pulse appears to be situated where a depression in an m = 2 mode should be, much like the averaged pulses in Figs. 14 and 15 sit within a depression. The frames (b)–(e) offer clear evidence that the intermittent Lorentzian pulses are generated by a nonlinear interaction between drift-Alfvén modes.

In the lower $I sat$ panels, the event appears as an eruption on the gradient of the filamentary bundle that travels a short distance around the structure in the direction of the drift waves before dissipating. As mentioned, this event is a reconstruction of an ensemble of events that occur at a similar location in the (x, y) plane as the reference probe; a reconstruction using the reference probe on the other side of the filaments produces a similar event. Thus, the occurrence of pulses of varying amplitude and width all around the filaments and also on the gradient at z₂ gives the impression that these eruption events occur sporadically all around the bundle and along the full axial length of the filaments. While the $I sat$ signal is a combination of density and temperature, there is strong evidence that the observed events correspond to anomalous energy transport. The blob moves into the cold background plasma somewhere between t = t₂ and t = t₃, by t = t₄ the $I sat$ signal is mostly gone and has dissipated within about 50 μs; presumably, this must be heat conduction along the field lines as this timescale is too fast for density outflows. In the future, conditionally averaged reconstructed sweeps^16,31 could be used to confirm this interpretation; additionally, the reference probes in this experiment were placed on the outer gradient to obtain information on the global drift-Alfvén modes, placing a reference probe in the middle of the bundle may allow for the negative events observed in the PDF at z₁ [Fig. 9(a)] to be spatially reconstructed in the same manner as presented here.

In this section, we relate the transport dynamics to deterministic chaos using a method known as the complexity–entropy plane or “CH-plane,” which can be used to distinguish between periodic, stochastic, and chaotic signals in a time series.³² Lorentzian-shaped pulses and exponential frequency spectra are signatures of chaotic dynamics and the CH-plane is another technique that can be used to confirm whether underlying chaotic dynamics are present. Maggs and Morales¹² made the first application of the CH-plane analysis to single filament experiments in the LAPD and demonstrated the chaotic nature of turbulent fluctuations on the gradient. The technique has since been applied to experiments involving magnetic flux ropes^33,34 and comparisons of laboratory turbulence and the solar wind.³⁵ Here we apply this technique to the full plane of Langmuir probe data in the multi-filament and single-filament cases.

The CH-plane is constructed from two measures of the signal, the permutation entropy and statistical complexity. Both measures may be calculated using a probability distribution introduced by Bandt and Pompe.³⁶ The distribution is calculated from the amplitudes of a signal of length N_t divided into $N t − D + 1$ sequential sets of length D, where D is referred to as the embedding dimension. The distribution represents the probabilities of the $N = D !$ permutations, ρ, of the ordering of the amplitudes of each set. For example, given the sequence $x = [ 1 , 4 , 5 , 1 , 3 , 2 ]$ and embedding dimension D = 3, there are N = 6 permutations and $N t − D + 1 = 4$ sets of length D: $[ 1 , 4 , 5 ] , [ 4 , 5 , 1 ] , [ 5 , 1 , 3 ]$⁠, and $[ 1 , 3 , 2 ]$⁠. These sets have the permutations ρ: (012), (120), (201), and (021), which each have probability $p ( ρ ) = 1 / 4$⁠; the remaining two permutations, (210) and (102), have no occurrences and zero probability. For a general time series ${ x t } t = 1 … N t$⁠, the probability for each ρ can be written,³⁶ 

where $#$ represents a count or “number of occurrences.” Once the set of N probabilities, $p ( ρ )$⁠, is known, the entropy and complexity are calculated to find the position of the signal in the CH-plane. The entropy is calculated using the Shannon entropy definition,³⁷ 

where S is the Shannon entropy, H_S is the normalized Shannon entropy, the sum is over all ρ, and $S max = S ( p e )$ where p_e refers to a uniform distribution across all ρ, i.e., $p e ( ρ ) = [ 1 / N , … , 1 / N ]$⁠. The statistical complexity is calculated using the Jensen-Shannon complexity defined as follows:³⁷ 

Central to applying the CH-plane technique to experimental data is the concept of sub-sampling or embedding delay.^12,35 The time series has a time step $Δ t$ determined by the sampling frequency, f_s. The CH-plane will resolve dynamical processes, be they stochastic, periodic, or deterministic, that occur on time scales greater than $2 Δ t$⁠, i.e., with a frequency below the Nyquist frequency, $f N = f s / 2$⁠. Real-world signals from dynamical systems are often characterized by noise beyond the frequencies of interest, such noise will dominate the probability distribution and suppress the ability of the CH-plane to extract evidence of chaotic dynamics. To handle the issue of noise, the signal is sub-sampled by an integer factor m to increase the time step between points in the permutation series ${ x t }$ to $m Δ t$⁠, effectively reducing the Nyquist frequency. Sub-sampling by taking every m sample will reduce the length of the time series to $N t ′ = N t / m$ and exclude the rest of the data; instead, here we use the length preserving technique³⁵ where the sub-sampling is accomplished by constructing sets of length D as ${ x t , x t + m , … x t + m ( D − 1 ) }$ with $t = 1 , . . N t − m ( D − 1 )$⁠, which preserves nearly all of the original data set for the calculation of $p ( ρ )$ and increases the sampling time.

In Fig. 17 is shown the CH-plane for three data sets: the 2D (x, y) planes at z₁ and z₂ investigated throughout this article in Figs. 17(a) and 17(b), respectively; and a 2D (x, y) plane of the single filament during the transport regime where it exhibits turbulent intermittent pulses.^10,12,19 Each data set has ten shots at each position with N_t = 8192, $Δ t = 1.28$ μs, length preserving sub-sampling with m = 4, and embedding dimension D = 5. The choice of $m Δ t = 5.12$ μs (Nyquist frequency $≈ 100$ kHz) on the order of τ for the Lorentzian pulses reduces the impact of noise while still resolving the intermittent events. First, a description of the different regions of the CH-plane is necessary. The range of the entropy of the signal is $0 ≤ H S ≤ 1$ but the complexity C_JS at a given entropy is bounded by a minimum and maximum (the black curves in Fig. 17); a description of the minimum and maximum conditions is given by Martin et al.³⁷ The left side of the plane is populated by periodic signals with low entropy. For a pure sine wave of frequency f with f_s approaching infinity $p ( ρ )$ would tend toward only two non-zero permutations, monotonically increasing and decreasing, and the entropy would approach $H S = ln ( 1 / 2 ) / ln ( N )$ (⁠ $≈ 0.145$ for D = 5, as in the figure). This is indicated by the blue curve that starts in the bottom left with signal to sampling ratio $f / f s = 1 × 10 − 4$ and goes up to $f / f s = 0.49$⁠. The lower right side of the plane is populated by stochastic signals with high entropy and low complexity. The dashed red line indicates fractional Brownian motion (fBm) with Hurst exponent,³⁸ H, ranging from $( 0 ≤ H ≤ 1 )$⁠, beginning in the bottom right. Signals below the fBm line are stochastic, and those lying above the fBm line in the upper part of the plane are characterized as chaotic signals. The locations of three chaotic maps, the Lorenz attractor,³⁹ Schuster map,⁴⁰ and Logistic map,⁴¹ are shown across the upper half of the plane; short descriptions of each chaotic signal and the parameters used to generate them are located in the  Appendix. Last, experiment shots of the afterglow plasma with no filament present and sampled the same as the experiment signals are shown in the bottom right corner; this demonstrates the stochastic nature of the probe signal in the absence of the filamentary structures.

In all three cases investigated in Fig. 17, there is clear evidence of chaotic dynamics; in the case of (c) this is expected from past results.¹² In Fig. 17(a), it is observed that the highest complexities occur at a radius between 0.9 and 1.3 cm, which corresponds to the peak of the drift-Alfvén wave modes and region of the largest intermittent events. Beyond r = 1.5 cm, the signals begin to approach the fBm line and stochastic processes are dominating. A similar outcome is seen at z₂ in Fig. 17(b) where the most complex signals are found between r = 0.2 cm and 0.7 cm. For the single filament (data collected at z₁) in Fig. 17(c), the chaotic dynamics are restricted to $r ≤ 0.5$ cm and for larger radii, the location in the CH-plane intersects the fBm line and indicates a strong stochastic nature. The radial extent of the chaotic dynamics in the tri-filament bundle is larger than for the single filament, even when taking into consideration the off-axis filament peaks. In the tri-filament case, the chaotic dynamics dominate up to $2 δ e$ from the profile peak while in the single filament case the dynamics are stochastically dominated beyond ∼1δ_e from the filament center. This result is in agreement with the enhanced extent of anomalous transport events in the interacting tri-filament bundle.

In this paper, the interrelationship between unstable turbulence and intermittent fluctuations in multiple interacting electron plasma pressure filaments has been investigated in a linear magnetized plasma column.

As an extension to our previous analysis of drift-Alfvén modes in the presence of multiple interacting filaments,²⁰ we have presented results on the axial variation of the plasma pressure profiles and drift-Alfvén modes in the filament bundle. An off-axis peak in the radial pressure profile occurs closer to the thermal sources with a shift of the pressure peak toward the center of the filaments at more distant axial locations. With this tendency, the maximum gradient of the pressure profile moves toward the filament center and similarly for the drift-Alfvén modes, which remain localized to the maximum pressure gradient region. Furthermore, through a mode decomposition procedure we have established that coupling between different azimuthal mode numbers occurs at both axial planes.

Following the radial and axial analysis of the drift-Alfvén modes, the probability distribution function (PDF) of the ion saturation current fluctuations (proportional to electron temperature and density fluctuations) was constructed, revealing their intermittent and non-Gaussian character. Near the center of the filament bundle the PDF follows a Gaussian distribution, with approximately equal numbers of positive and negative fluctuations about the mean. At the off-axis peak in the radial pressure profile, there is a weak enhancement of the negative fluctuations about the mean indicating a tendency for depletions of temperature and density. However, in the outer gradient region, the fluctuations that predominate are positive compared to the mean. The largest positive skewness of the PDFs coincides with the peak in the drift-Alfvén mode amplitudes.

In the context of this experiment, intermittent fluctuations are linked to the presence of Lorentzian pulses. The presence of Lorentzian pulses in the underlying time series leads to an exponential frequency power spectrum which is a signature of deterministic chaos. Using a large collection of Lorentzian pulses in the time series, we have made a statistical characterization, both temporally and spatially. Based on the PDF of the amplitude of the intermittent pulses at different radii and two different axial locations, it is found that the pulse amplitudes decrease away from the plasma pressure gradient. For the axial location closer to the source, the pulse width increased with decreasing pulse amplitudes, whereas at the axial location further from the source where the filaments are merged, the pulse width was nearly constant for increased pulse amplitude. By using a cross-conditional averaging method involving a nearby reference probe, it has been possible to make a spatial reconstruction of the pulse event. This establishes the spatial scale of the pulse, roughly twice the electron collisionless skin depth in both the radial and azimuthal directions, and associates it with the underlying drift-Alfvén mode pattern. The intermittent pulses are mainly blobs of temperature moving radially and azimuthally away from the filaments and are the dominant source of cross field energy transport.

Another approach to infer whether a physical system contains deterministic chaos in the underlying dynamics is to construct the complexity–entropy (CH-plane) from the fluctuation time series. Following the standard procedure of permuting the fluctuation amplitudes to form the permutation entropy and computing the statistical complexity, we have determined the radial extent of the chaotic dynamics region and demonstrated that they are mainly localized to the electron plasma pressure gradients. These coincide with the source region of the intermittent pulse events and the peak in the amplitude of the saturated drift-Alfvén modes.

In conclusion, we have presented observations of strongly intermittent turbulence in multiple interacting plasma pressure filaments. This system has underlying deterministic chaotic dynamics that exhibit intermittency through Lorentzian pulses associated with coherent structures evolving from drift-Alfvén mode instabilities. It has been demonstrated that large amplitude temperature fluctuation bursts are related to spatiotemporal coherent turbulent structures that propagate radially outward as well as azimuthally, away from the filamentary structures. Therefore, the cross field energy transport is not simply characterized by a quasi-Gaussian diffusive model but better described by chaotic and intermittent transport events, which contain the physical effects of electron pressure gradient-driven drift-Alfvén instabilities in the presence of nonlinear convection of density, temperature, and vorticity.

The authors would like to thank Professor Troy Carter for helpful discussions about the manuscript. Authors R.D.S. and T.S.G. acknowledge support from the Natural Sciences and Engineering Research Council of Canada (NSERC). Author S.K. acknowledges support from NSERC, NASA Grant Nos. 80NNSC19K0264 and 80NNSC19K0848, and NSF Grant No. 1914670. The experiments were performed at the Basic Plasma Science Facility supported by DOE and NSF, with major facility instrumentation developed via an NSF Award No. AGS-9724366.

The authors have no conflicts to disclose.

Scott Karbashewski: Conceptualization (equal); Data curation (lead); Formal analysis (lead); Investigation (equal); Methodology (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (equal). Richard Sydora: Conceptualization (equal); Formal analysis (supporting); Funding acquisition (lead); Investigation (equal); Methodology (supporting); Project administration (lead); Writing – original draft (supporting); Writing – review & editing (equal). Bart Van Compernolle: Conceptualization (equal); Data curation (supporting); Formal analysis (supporting); Investigation (supporting); Resources (lead); Writing – review & editing (supporting). Thomas Simala-Grant: Data curation (supporting); Formal analysis (supporting); Writing – review & editing (supporting). Matthew Joseph Poulos: Conceptualization (equal); Investigation (equal); Writing – review & editing (supporting).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

In this appendix, we outline the chaotic signals and associated parameters used in the analysis of the complexity–entropy plane. These are all common chaotic systems with many resources available for further reading, and only a basic description of each is given here.

The Lorenz System is a system of ordinary differential equations with solutions that exhibit chaotic behavior for certain values of the constant parameters within the system; for parameters giving these chaotic solutions, it is referred to as the Lorenz Attractor. The system is given by the equations,³⁹ 

where the system parameters ρ = 40, σ = 10, and $β = 8 / 3$ are used in our analysis to yield chaotic solutions.

The Logistic map is a simple discrete recurrence relation that exhibits chaotic behavior for certain parameters and is defined as⁴¹ 

and the system parameter r = 3.9076 is used to generate the chaotic signal used in our analysis.

The Schuster map is another discrete map that exhibits chaotic behavior described by the relation,⁴⁰ 

where the parameter Z = 3/2 yields the chaotic series used in our analysis.