Observation of Kolmogorov turbulence due to multiscale vortices in dusty plasma experiments

Laboratory dusty (complex) plasmas comprise of electrons, ions, neutral gas particles, and highly charged dust particles.^1,2 Dust particles generally have sizes up to tens of $μ$m and typically attain a high negative charge of $104$ e due to the high mobility of electrons. However, under atypical conditions, the dust particles can also acquire positive charges transiently or stably, as reported in various afterglow dusty plasma experiments.^3,4 Their dynamics are slow (a few tens of Hz) compared to electron and ion timescales. The length scales usually range from a few millimeters to a few centimeters.¹ The slow dynamics and convenient length scales render it an ideal experimental platform for accurately monitoring complex physics phenomena. Some examples are phase transition, longitudinal and transverse waves, solitons, shocks, vortices, and complex flows.^5–13 

Dusty plasma experiments with precision imaging diagnostics enable the study of complex flows. It includes two- and three-dimensional laminar, sheared, rotating, chaotic, and turbulent flows.^7,9,14–16 These low Reynolds number flows can be classified as either hydrodynamic or visco-elastic. This classification depends on the interaction strength of the charged dust particles relative to their thermal energy in the plasma. Dusty plasma flows are often observed in experiments due to the interplay of various forces acting on dust particles. These flows can be tracked from particle scales to hydrodynamic scales.^17–19 In this study, we observed the turbulent flow caused by the dust vortex, particularly within the diffuse region of a DC glow discharge plasma. We found the alignment with Kolmogorov power law scaling in both spatial and temporal domains of the energy spectrum, thus demonstrating Kolmogorov turbulence in afterglow-like dusty plasma conditions.

Dusty plasma experiments provide a unique platform for studying turbulent flows across both macroscopic (system size) and microscopic (average interparticle separation) scales. This precision is achieved by tracking individual dust particles' trajectories, even in the liquid phase. These experiments serve as an alternative to Navier–Stokes equation-based numerical models, which depend on the continuum approximation. Such models are computationally expensive because they require numerous grid points to resolve scales close to the particle level.^20,21

Another notable aspect of dusty plasma experiments is that they can easily be formed in a strongly coupled regime. Turbulent mixing within such a regime holds significance in astrophysical flows, inertial confinement fusion, and other areas.^22–26 Furthermore, the strong coupling phenomenon is recognized for introducing viscoelastic features in dust flows.²⁷ Consequently, these experiments serve as a foundation for investigating turbulence within the low Reynolds number regime, often called viscoelastic turbulence.²⁸ Finally, dusty plasma experiments provide a test bed to observe the effect of gravity on turbulent mixing. Theoretically, gravity is predicted to steepen energy transfer scaling in at least two-dimensional flows.²⁹ In our experiment, the dust flows, and rotation forms naturally at a particular set of discharge conditions. We will examine the observed flow for turbulence characteristics based on the possibilities mentioned earlier.

Nonlinear mixing and turbulence in dusty plasma can broadly have two origins: nonlinear wave interaction and vortex formation.^5,30–36 In numerical modeling, turbulent mixing in dust fluid occurs by initiating energy/perturbation in vortices of system size or by introducing various seed instabilities.³⁷ The main question is whether the turbulence characteristics in strongly coupled plasmas resemble or differ from those in hydrodynamic fluids. In the context of turbulence in hydrodynamic fluids, Kolmogorov proposed that the energy spectrum $E(k)$ follows a $k−5/3$ scaling in the inertial range. It applies to any homogeneous, isotropic turbulent flow forced at a large scale. Here, k is the wavenumber associated with the particular eddy in the turbulent fluid system, and $E(k)$ is the energy associated with the respective eddy.³⁸ 

For dusty plasmas, Schwabe et al.³⁰ experimentally demonstrated a dual cascade in the spectrum. They found that $E(k)∝k−5/3$ for $k<kl$ and $E(k)∝k−3$ for $k>kl$⁠, where $kl$ is an intermediate forces wavenumber. Their results are consistent with theoretical predictions of two-dimensional forced turbulence.^38–40 Tsai and Lin-I⁴¹ demonstrated that the dust acoustic wave turbulence exhibits $f−2.6$ scaling in a frequency spectrum. Using kinetic numerical simulations, Joshi et al.³⁴ demonstrated turbulent mixing due to the vortex motion in the dust flow past the obstacle. Their results follow the hydrodynamic Kolmogorov turbulence scaling of $−5/3$ in their energy spectrum in both the time and spatial domains. Using molecular dynamics simulation, Wani et al.⁷⁸ demonstrated the turbulence in a dusty plasma system using Rayleigh–Taylor instability as an initial condition. Their spectra align with the Bolgiano–Obukhov theory, a typical two-dimensional Rayleigh–Taylor turbulence. While all previous experimental works report wave turbulence, this work is the first to demonstrate fully developed Kolmogorov turbulence in dusty plasma due to multi-scale vortex flow. The vortex flow appears stationary at the largest scale of energy infusion. We could show the scaling of the energy cascade close to $−5/3$ in the wavevector and frequency domains.

In this paper, we report a few large-scale vortices and small-scale eddies in a dust cloud. It levitates in the diffused plasma far from the electrode in a DC discharge. The dust cloud remains stable for several hours due to the balance of electrostatic $(qE→)$ and gravitational $(mg→)$ forces. This stability is maintained as long as an external source continuously provides energy, and the dust cloud resides in an energy bath of background plasma. The dust cloud is 3D in shape. However, the image analysis and particle image velocimetry (PIV) are carried out on the basis of the two-dimensional imaging of the cloud. We estimate two-order-of-magnitude background plasma density variations between the electrode and the diffused plasma region. The observed vortex is reproducible on the plasma device, Shivalik Plasma Device-I (SPD-I) device.⁴² We presented the turbulence features in the rotating dust cloud of the diffused plasma region. Understanding the flow dynamics in diffused plasma regions is interesting due to their relevance to the interstellar universe and the lunar surface, which have minor density variations. We observe a scaling consistent with Kolmogorov's $−5/3$ in the frequency and spatial domains of the energy spectra. Also, in all energy spectra calculations, the kinetic energy is computed using only the fluid velocity components without incorporating the masses of fluid elements. To further support our results, we have demonstrated the second-order structure functions following the Kolmogorov law and the probability distribution function (PDF) of velocity gradients.

The paper is organized as follows: Sec. II shows details of the device, image processing technique, and estimated plasma parameters. Section III provides observations and analysis, followed by an associated discussion. Finally, Sec. IV summarizes the outcomes with future possibilities.

The experiments are performed on a newly commissioned laboratory plasma device, Shivalik Plasma Device-I (SPD-I).⁴² This plasma chamber is a horizontally laid (Z-axis) glass cylinder with two asymmetric electrodes kept vertically in a parallel plate arrangement, as shown in Fig. 1. The glass chamber is electrically floating using Teflon connectors, with the top electrode as an anode and the bottom grounded electrode as a cathode. In the present experiments, we separated both electrodes by 7.5 cm. Figure 1 depicts this DC discharge-based experimental setup. The details on configuration, including the device dimensions, electrodes, and diagnostics, are in Sharma et al.⁴² We introduced argon, the working gas, into the vacuum chamber through a port adjacent to the pumping port to prevent any directed flow of neutral gas into the experimental region of interest. This approach ensures that only diffused gas enters the system.⁴³ We repeatedly purged the device with argon gas at high pressures (⁠ $∼5$ mbar) for at least 10 minutes before maintaining working pressure to conduct the experiment. The plasma was formed in an argon gas background at a low working pressure of 0.12 mbar.

Polydispersive kaolin particles, ranging in size from 4 to 8  $μ$m, are dispersed over the cathode before creating the vacuum. We performed a slow initial pumping to ensure the Kaolin particles remained inside the device. As soon as plasma creation takes place, the dust particles covering the cathode, along with some of its floating fraction, become charged. A void quickly forms over the cathode because of the competition between the ion drag force and the radial electric field. During the void formation, the electrode expels almost all floating dust particles. However, these dust particles then form a 3D cloud in the diffused plasma region over the extended sheath on the glass surface. The 3D cloud forms, possibly for two reasons: first, because of polydispersive Kaolin particles, and second, because of the 3D potential well formation in the diffused plasma region. The second reason is the predominant cause of the 3D dust cloud in the present experiment. We validated it using monodispersive MF particles in the SPD-I device⁴² and a similar 3D dust cloud formed. Three forces may contribute to this potential well formation: (1) A radially outward force on dust due to void formation over the cathode; (2) the force generated by the lower glass surface's sheath electric field; and (3) the force emanating from the electric field profiles of the end-flange sheath.

A green LASER with a 532 nm wavelength and a variable power capacity of up to 300 mW illuminates the dust particles. A cylindrical lens is placed before the LASER to convert it into a thin sheet of thickness 150  $μ$m. A high-speed camera captures the illuminated dust particles at a 90-degree angle to the LASER sheet (X-Z plane). We have also made a few observations by keeping the LASER sheet in the Y-Z plane to support the three-dimensional nature of the dynamics. We have used two different cameras. The PCO Edge 5.5 is utilized for recording up to 100 fps with a resolution of 2560  $×$ 2180 pixels, chosen for its superior spatial resolution. The Phantom v2512 records up to 1000 fps with a resolution of 1064  $×$ 1064 pixels, selected to achieve a superior temporal data profile. PIV has served as the primary tool for extracting velocity fields from the recorded image sequences (detailed in Subsection II B). We also extract particle positions and velocities by tracking particles over time for at least two experimental datasets. We used these particle datasets to validate PIV results and predict the static structural properties of dust clouds.

We conducted multiple campaigns to study flows and turbulence to ensure the reproducibility and consistency of observed dust rotation phenomena. A total of six independent experimental campaigns were conducted at the same working pressure of 0.12 mbar, as listed in Table I. All of these campaigns have been performed using polydispersive kaolin microparticles. The schematic in Fig. 2 shows tentative dust rotation locations from different campaigns. These campaigns successfully captured the X-Z plane of the dust cloud, ensuring that the laser sheet remained at the center of the glass flange at location $y=0$⁠. The hemispherical, shaded region represents the dust void. In each campaign, we observed dust rotation with variations in flow characteristics. Specifically, I, II, III, and VI form in reasonably different plasma environments. We anticipate an interplay of different sets of forces^44–48 to be responsible for individual vortices in the extended diffused region. It is interesting to note that in all these conditions, the nature of turbulence almost follows the standard forced Kolmogorov theory.

Table II enlists parameters for the background plasmas and the dust cloud. The electron temperature, $Te$⁠, and the plasma density, $ne,i$ (in quasineutral plasma), were measured using the single Langmuir probe above the electrode and in the diffused plasma region in the absence of dust.⁵² Over the cathode, $Te$ is in the range of $1.5–4$ eV from the center to the edge. In the diffused region where dusty plasma typically forms, $Te$ is relatively high, in the range of 3 to 7 eV at different locations within the domain span. The plasma density varies by up to two orders of magnitude from the electrode to the diffused region of interest, as shown in the table.

We used a microscope to obtain the kaolin particles' size distribution. Most particles ranged from 5 to 8  $μ$m, although there were a few outliers with sizes as large as 50  $μ$m. Based on the known mass density of kaolin used in the experiments, the typical mass of dust particles ranged from $1.6×10−13$ to $6.7×10−13$ kg. Out of the six campaigns, we could track particle positions with reasonable accuracy using the ImageJ software⁵¹ for two campaigns. We used the particle data from Fig. 3(b) for 145 successive frames and a particle tracking algorithm to gather particle velocities (in mm/s). Determining the accurate average interparticle separation, $(a)$⁠, for dust clouds composed of particles of varying sizes and experiencing different spatial plasma environments is challenging. However, we calculated “a” from regions near vortex motion with visually no flow but only thermal motion. We obtained an estimated value of $a≈224 μ$m and a typical dust density, $nd≈2.3×104$ cm⁻³, in the diffused plasma region. Figure 3(b) depicts one such particle arrangement (from campaign V) we used to estimate a.

Dust particles typically have a charge in the range of [2–5] $×104e$⁠. The floating potential on a spherical dust particle is the basis for the collisionless orbital motion limited model, which is used to compute the charge. We fitted the theoretical Gaussian velocity distribution to the probability distribution function (PDF) of the particle velocity data (averaged PDF from 145 frames), which we plot in Fig. 3(a). A theoretical fitting in comparison with experimental data yields dust temperature ranges of 0.6–2.2 eV. We see that in a three-dimensional situation, the dust temperature will vary since it is computed using a two-dimensional velocity field. For two datasets (II and VI), we additionally computed the 2D radial distribution function (RDF). The dynamics are thermal, and there is no bulk flow in the area used for this computation. With a tiny initial peak, the RDF in Fig. 4 illustrates the equilibrium dust dynamics in a somewhat liquid regime.

Image processing and PIV are crucial for accurately extracting fluid velocity in the dust cloud. We detail the processes adopted from recorded videos to ascertain velocities over a square grid. We processed video data from six campaigns (I–VI) to study fluid velocities and turbulent mixing. Campaigns I and II were recorded at 25 fps, III and IV at 101 fps, and V and VI at 125 fps. All original images were pre-processed using the open-source ImageJ software by converting them into grayscale. Thresholding was then applied to convert all pixels with an intensity range of 50–220 to maximum brightness, while other pixel intensities were set to zero. A convolution filter sharpened the particles for better clarity, especially given the significant particle-size asymmetry. These data were then ready to export for PIV studies.

We used MATLAB-based PIVlab software for the PIV analysis of the pre-processed image data.⁵³ After importing the images, we selected the square region of interest (RoI) to study vortex and turbulence flows. Unwanted regions (blank spaces with no dust) were removed by masking. For the PIV analysis, we used a $32×32$ pixel interrogation area with a 75% overlap of windows. The Fast Fourier Transform (FFT)-based cross correlation algorithm calculated velocities (pixels/frame) for all windows within the RoI. The velocity data were then calibrated based on actual experimental dimensions, using pre-marked references in the images for length and fps for timescales. Unwanted high-velocity data significantly higher than the bulk flow speed were removed using velocity histogram cropping. Finally, we obtained N-1 sets of velocity field data from N images and used them to profile velocity fields, vorticity, and energy spectra.

We report the observation of self-excited dust rotation dynamics in the 3D dust cloud within a diffused plasma. We analyze its properties to understand its turbulence characteristics. The following is a detailed discussion of the results.

The dust cloud in the diffused plasma region showed a variety of collective dynamics as we varied the discharge voltage between 230 and 600 volts. As listed in Table I, the self-sustained dust rotation has been observed for values at both edges of the voltage range. Figures 5 and 6 show time-elapsed images for different campaigns where dust rotation is visible from particle trajectories. The inset in Fig. 5(1) shows the tracking of selected particles in the box region. The particle tracking was performed using the ImageJ software's tracking algorithm. The exact location of each dust cloud with reference to the center of the cathode (in the Z-direction) and the surface of the glass (in the X-direction) is mentioned as axis labels. All other dust rotations are observed typically at the same location and similar discharge parameters except for campaigns I and II. For each campaign, the vortex motion is quite stable, and dust continues to rotate until we change the plasma creation conditions. We have provided dust rotation movies for each campaign as the supplementary material.⁵⁴ 

We understand that a typical dust rotation in an unmagnetized dusty plasma may arise due to counteracting forces at spatially different locations.¹³ While gravity is always one force acting downward, the radially outward ion drag force in the void, the spatially varying forces due to the sheath electric field are few possibilities in the counter direction, or at least their one component is opposite to the gravity.

A comprehensive analysis to establish the cause of dust rotation is beyond the scope of this paper. However, we ruled out some of the possibilities. A temperature gradient in neutral gas can give rise to thermophoretic force, which can affect dust dynamics.² Typically, in glow discharge, the temperature of the cathode can increase reasonably due to ion bombardment after a few hours of operation.⁴⁴ In our experiment, we observed the dust vortex soon after reaching the operating discharge voltage range before the cathode was reasonably heated. This suggests that thermophoretic force may not influence the vortex dynamics significantly. We also ruled out the possibility of a neutral flow in the vacuum chamber affecting the dust dynamics. We repeated the vortex experiment by simultaneously closing the pump and the Ar gas inlet to minimize any kind of neutral flow in the vacuum chamber. In this condition, too, we observed the dust vortex, thus ruling out any role of neutral flow.

In earlier dusty plasma experiments with vortex observation, the causes involved variations in microparticle sizes,^55,56 charge gradients,^57,58 or a non-zero curl of the forces exerted by the plasma on the microparticles.^59–61 In the future, we plan to study the force balance comprehensively in our system. The goal of this paper is limited to demonstrating turbulence in a dusty plasma medium.

Figures 7 and 8 provide flow velocity vectors and vorticity plots for all campaigns. The vorticity plots were obtained after taking the mean of velocity field profiles from 100 frames, which profoundly helped to visualize the large-scale rotations. The averaging cleans up the fluctuations and reflects the large-scale vortex structures. However, these fluctuations are important when evaluating the energy spectrum in Sec. III B. The vorticity profiles of all campaigns reflect that the dust rotation is clockwise for campaigns III–VI as shown in Fig. 8, and the same is counterclockwise for campaigns I and II (see Fig. 7). This is observable from the velocity flow vectors and the colors of the vorticity profiles in different plots.

An important feature of turbulent flows is that they exhibit multi-scale energy cascades, which are typically demonstrated through the spatial energy spectrum.^62,63 Figure 9 presents subplots displaying the energy spectra for each campaign, which demonstrate the distribution of kinetic energy across different scales or wavenumbers. The spectra are calculated in the following way: For each recorded dataset, we take the PIV and get multiple velocity-flow profiles for each campaign. Once we have a coarse-grained flow velocity $u(x,t)$ for a given snapshot, we perform a Fourier transform on $u(x,t)$ to obtain $u(k,t)$⁠, from which we calculate the modal energy $E(k,t)=|u(k,t)|2/2$⁠. After that, we find the one-dimensional shell spectrum $E(k,t)$ by adding up the modal energies within a certain range of $k−1≤k′≤k$⁠, as described by Frisch and Kolmogorov.⁶³ 

We estimated the energy spectra for each flow profile of the recorded datasets separately. Then, we plotted the average energy spectra over several time frames with error bars, as seen in Fig. 9 for all the campaigns. The wavevectors are normalized with $kmax$⁠, where $kmax$ is the maximum wavenumber associated with the system. The values of $kmax$ for each campaign are 6.334, 5.69, 5.69, 5.69, 4.91, and $4.91mm−1$⁠, respectively. Note that the units of $kmax$ are $mm−1$⁠, and the units of the energy spectra are $(mm/s)2$⁠. We perform normalization to present the spectra on the same scale without losing generality.

We observe that the energy fed at the large scales cascades to the smaller scales following the standard -5/3 scaling for all campaigns, which is consistent with Kolmogorov's theory.⁶⁴ The corresponding scaling in the real space is observed by analyzing the second-order structure function discussed in Sec. III C. However, note that the wavevector ranges vary slightly for different campaigns. Also, the total energy range in the spectra varies up to an order of magnitude as we look for the spectra for I–VI. It is based on the vortex domain size and finite flow velocities. For example, in Fig. 7, the vortex size of I is smaller than that of II; hence, its energy spectrum has a larger magnitude range.

The flow dynamics is three-dimensional, but the present analysis is carried out only in the x-z plane due to the limitation of recording only one 2D plane using a laser sheet and a camera. However, we have individual movies of the y-z and x-z planes (asynchronous) for campaigns V and VI. Here, we provide in Fig. 14 (see the  Appendix), the energy spectrum and the vorticity plots for campaign V in the y-z plane. These plots demonstrate the 3-D characteristics of turbulent vortex flow in the cloud at the specific position of $x=5$ mm. The least-square fitted scaling (the slope) obtained is $−1.595±0.099$⁠, consistent with our claims of Kolmogorov scaling within experimental uncertainties.

We analyze the structure function for examining the statistical characteristics of turbulent flows.⁶³ The structure functions display the moments of velocity increase. $δurp=⟨[u(r→+l)−u(r→)]p⟩$⁠, where r is a reference point, p is the order of the structure function, and l is the distance between two points in a velocity field.⁶⁵ Based on Kolmogorov's theory (Kolmogorov 1941 b), for fully developed homogeneous turbulence, the scaling relationships in the inertial range of turbulence should follow the form $δurp∝rξ(p)$⁠, where $ξ(p)=p/3$ roughly.⁶⁶ 

In this study, we calculate the second-order structure function for all datasets by determining the difference in velocity between a reference point r and varying separation distances l, as shown in Fig. 10. To maintain generality, we computed the structural functions using a time average of 90 velocity frames. The fitted slopes yield the scaling exponent close to 0.66 with standard deviation, which is almost consistent with Kolmogorov theory.³⁸ We also see that the $S2(l)$ values for campaigns III and IV are two orders of magnitude lower. This is because the velocity order is lower in these vortex regimes, which is also shown in the energy spectra, Fig. 9.

Using analytical arguments, we also estimated the Kolmogorov dissipation length and time scales. These scales were calculated based on the mean energy transfer rate $ϵ∼vl3/l$⁠, where $vl=2 mm/s$ is the average velocity and $l=a$ is the average interparticle separation. The local viscosity is $ν∼vll$⁠. With $ϵ$ and $ν$⁠, we computed the Kolmogorov length scale as $ld=(ν3/ϵ)1/4=110 μm$ and the velocity scale as $vd=(νϵ)1/4=2 mm/s$⁠. Consequently, the dissipative timescale is $τd=ld/vd=0.0478 s$⁠, leading to a dissipation frequency of approximately $21 Hz$⁠. This is also clear from Fig. 11, which shows the dissipation frequency at almost $21 Hz$⁠.

Since the turbulent velocity fields exhibit significant fluctuations, it is essential to comprehend their statistical properties. Consequently, researchers have focused on studying velocity fields' PDFs and their derivatives. Velocity PDFs provide insights into the distribution of flow velocities at different points in space.^67,68 The PDF of velocity fields should follow a Gaussian distribution in fully developed turbulent flows, but their derivatives often deviate from the Gaussian distribution, especially at the tails.⁶⁹ 

Figure 12 presents the PDFs of velocity gradients, $δvz$ and $δvx$⁠, which deviate from the Gaussian distribution (shown by the black dashed lines) at the tails, indicating intermittency. Note that the PDFs shown are averaged over 90 datasets for each campaign. The fitted plots are obtained by providing $δvx$ and $δvz$ dataset to the standard MATLAB function “fitdist” for “normal” distribution. The distribution returns a mean value and the associated standard deviation in each case. Table S2 detailing mean velocity, standard deviation, skewness, and kurtosis is included as the supplementary material.

In this work, we demonstrated the presence of Kolmogorov turbulence in a three-dimensional dusty plasma. The dust cloud formed naturally as a result of plasma and device configuration forces. The energy spectrum scaling of −5/3 and the second-order structure factor with a slope of 2/3 are reproducible for different campaigns, regardless of the vortices' origin and location. We could also demonstrate the characteristics of stationary turbulence due to −5/3 scaling in both spatial and time-series spectra. The self-consistent nature of the vortex motion of charged dust clouds also points to constant energy pumping up at the largest scales and dissipating as a thermal motion of dust particles. Overall, the dusty plasma medium acts as a driven dissipative system.

An important aspect of the present work is the demonstration of fully developed turbulence. Traditional understanding suggests the possibility of fully developed turbulence at a high Reynolds number $Re$⁠.⁷⁰ In dusty plasmas, earlier studies reported wave turbulence at low $Re$⁠.^71,72 The present work examines and reports fully developed turbulence through multi-scale vortex dynamics. Our calculations suggest typical values of $Re$ between 8 and 50. Section V appendix includes the calculation of $Re$⁠. These results suggest dusty plasmas' similarity to viscoelastic fluids where elastic turbulence is observed at small $Re$ values. This anticipation contradicts presented scalings from energy spectra and structure functions that suggest close alignment with Kolmogorov's prediction. This remains an open question, and we will try to examine this in future studies.

See the supplementary material for the following: (I) we provide an agreeable picture of full three-dimensional vortex flow using three different x-z planes at locations y=−1,0,1 mm; the energy spectra from all three planes in Fig. S1; (II) the movies for campaigns I–VI provided in '.m4v' format; and (III) the statistical parameters associated with the velocity gradients for campaigns I–VI.

The authors acknowledge Professor Mahendra Verma for his valuable discussions. ST acknowledges support from IIT Jammu PRAISE grant PRA-100003 for the experimental setup. ST and RW acknowledge partial support for this work through SERB Grant No. CRG/2020/003653. SB was supported by the U.S. Department of Energy under Contract No. DE-AC02-09CH1146.

The authors have no conflicts to disclose.

Sachin Sharma: Conceptualization (equal); Data curation (lead); Formal analysis (equal); Investigation (equal); Methodology (equal); Software (lead); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). Rauoof Wani: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Supervision (equal); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). Prabhakar Srivastav: Data curation (supporting); Validation (supporting); Visualization (supporting); Writing – review & editing (supporting). Meenakshee Sharma: Data curation (supporting); Validation (supporting); Visualization (supporting); Writing – review & editing (supporting). Sayak Bose: Conceptualization (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). Yogesh C. Saxena: Validation (supporting); Visualization (supporting); Writing – original draft (supporting); Writing – review & editing (supporting). Sanat Tiwari: Conceptualization (lead); Funding acquisition (lead); Investigation (equal); Methodology (equal); Project administration (lead); Resources (lead); Supervision (lead); Validation (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal).

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

We calculated the Reynolds number, $Re$⁠, for our experiments using the following parameters: $Γ=1500$⁠, $κ=1.5$⁠, $m=6×10−13$ kg, $n=2×1010$ m⁻³, $ρ=1.2×10−2$⁠, $a=250μ$m, and $ωpd=60$ Hz. From the literature, we took the normalized dynamic shear viscosity for 3D dusty plasma to be 0.2.^73–75 

The normalization factor is $η*=mnωpda2=4.5×10−8$ Pa.s. Thus, the dynamic shear viscosity is $η=0.2×η*=9×10−9$ Pa.s. Using a typical length scale $l=2×10−2$ m and a typical velocity scale $v=1×10−3$ m/s, we calculated the Reynolds number as $Re=ρlv/η$ = 26. Depending on the range of parameter values, the Reynolds number in these experiments varies between 8 and 50.

Because of intermittency,⁷⁶ the Kolmogorov $p/3$ scale in structure functions deviates at higher orders. In particular, the intermittency in turbulent flows is captured by the sixth-order structure function. According to Refs. 76 and 77, the scaling exponent for the sixth-order structure functions differs from Kolmogorov's theory and is $(p/3−μ)$⁠. The sixth-order structure functions are plotted in our work, as Fig. 13 illustrates. The intermittency parameters for campaigns I, II, and V are determined to be $μ=0.73$⁠, 0.86, and 1.2, correspondingly. Nevertheless, the parameters for campaigns III, IV, and VI are $−0.5$⁠, $−1$⁠, and $−0.49$⁠; for these campaigns, we need to perform more research as we do not yet have a satisfactory explanation.

We also provide the energy spectrum plot for campaign V in the Y-Z plane. The Kolmogorov scaling is approximately visible in the intermediate range of wavenumbers (see Fig. 14). Note that the approximate location of the dust cloud is at x = 5 mm from the Y-Z plane.