Unified transport scaling laws for plasma blobs and depletions

TL;DR

This paper derives and validates unified scaling laws for plasma blob and depletion dynamics, revealing how size, amplitude, and compressibility influence their velocities and accelerations in a gravitational field.

Contribution

It introduces a comprehensive empirical model that unifies the scaling laws for plasma blobs and depletions, accounting for compressibility effects and matching numerical simulations.

Findings

Radial velocity scales with the square root of size and amplitude for incompressible flows.

Compressibility alters the velocity scaling, making it linear with amplitude for certain ratios.

Depletions accelerate faster than blobs due to reduced inertia.

Abstract

We study the dynamics of seeded plasma blobs and depletions in an (effective) gravitational field. For incompressible flows the radial center of mass velocity of blobs and depletions is proportional to the square root of their initial cross-field size and amplitude. If the flows are compressible, this scaling holds only for ratios of amplitude to size larger than a critical value. Otherwise, the maximum blob and depletion velocity depends linearly on the initial amplitude and is independent of size. In both cases the acceleration of blobs and depletions depends on their initial amplitude relative to the background plasma density, is proportional to gravity and independent of their cross-field size. Due to their reduced inertia plasma depletions accelerate more quickly than the corresponding blobs. These scaling laws are derived from the invariants of the governing drift-fluid equations and agree excellently with numerical simulations over five orders of magnitude. We suggest an empirical model that unifies and correctly captures the radial acceleration and maximum velocities of both blobs and depletions.

Equations

$$\frac{\partial n}{\partial t}+\nabla \cdot \left(n{\mathbf{u}}_E\right)$$

$$\nabla \cdot \left(\frac{n}{\Omega }\left(\frac{\partial }{\partial t}+{\mathbf{u}}_E\cdot \nabla \right)\frac{{\nabla }_{\perp }\varphi }{B}+n{\mathbf{u}}_d-n{\mathbf{u}}_g\right)$$

$$\frac{\partial \rho }{\partial t}+\nabla \cdot \left(\rho {\mathbf{u}}_E\right)+\nabla \cdot \left(n({\mathbf{u}}_{\psi }+{\mathbf{u}}_d+{\mathbf{u}}_g)\right)$$

$$\frac{\mathrm{d}}{\mathrm{d}t}\left[T_{e}S(t)+H(t)\right]=0,$$

$$\frac{\mathrm{d}}{\mathrm{d}t}\left[E(t)-G(t)-H(t)\right]=0,$$

$$(MV{)}^2$$

$$\le 2{\left(\int \mathrm{dA}{\left[n\ln (n/n_0)-(n-n_0)\right]}^{1/2}\sqrt{n}{\varphi }_y/B\right)}^2$$

$$\le 4S(0)E(t)/m_i$$

$$n(\mathbf{x},0)$$

$$\frac{A_0}{g}=\mathcal{Q}\frac{2S(0)}{M}\approx \frac{\mathcal{Q}}{2}\frac{△n}{n_0+2△n/9}.$$

$$(MV{)}^2\le 4T_{e}S(0{)}^2/m_i.$$

$$\frac{\max |V|}{C_{\mathrm{s}}}=\mathcal{Q}\frac{2S(0)}{M}\approx \frac{\mathcal{Q}}{2}\frac{|△n|}{n_0+2/9△n}\approx \frac{\mathcal{Q}}{2}\frac{|△n|}{n_0}.$$

$$\frac{\max |V|}{C_{\mathrm{s}}}=\mathcal{R}{\left(\frac{\ell }{R_0}\frac{|△n|}{n_0}\right)}^{1/2}.$$

$$M_{\text{i}}\frac{dV}{dt}=(M_{\text{g}}-M_{\text{p}})g-c_{1}V-\mathrm{sgn}(V)\frac{1}{2}{c}_2{V}^2.$$

$$c_1=\pi {\ell }^2{n}_{0}g/C_{\mathrm{s}}.$$

$$c_2=\mathcal{Q}\pi n_0\ell /{\mathcal{R}}^2.$$

$$\frac{\max |V|}{C_{\mathrm{s}}}=\left(\frac{{\mathcal{R}}^2}{\mathcal{Q}}\right)\frac{\ell }{R_0}\left({\left(1+{\left(\frac{\mathcal{Q}}{\mathcal{R}}\right)}^2\frac{|△n|/n_0}{\ell /R_0}\right)}^{1/2}-1\right)$$

Full text

Fluctuation induced transport across magnetic field lines is ubiquitous in magnetized plasmas in various conditions. In the scrape-off layer of tokamaks field aligned plasma pressure perturbations universally appear. These perturbations are spatially localized when viewed in a plane perpendicular to the magnetic field and are often referred to as blobs. They mediate a significant amount of the radial particle and energy flux on plasma facing components and thus critically determine their lifetime Krasheninnikov (2001); Antar et al. (2001); D’Ippolito et al. (2002); Boedo et al. (2003); Garcia et al. (2005, 2006); Myra et al. (2006); Theiler et al. (2009); Carralero et al. (2015). Recent efforts in stochastic modeling relate the radial density profiles of magnetically confined plasmas to the amplitude, size and radial velocity of individual uncorrelated transport events such as blobs Garcia et al. (2016). Analysis of experimental data support the predictions of this stochastic model: probability density functions, auto correlation and power spectra as well as threshold level crossings of the turbulent fields are in good agreement with theoretical predictions Garcia (2012); Garcia et al. (2013, 2015); Theodorsen et al. (2016); Kube et al. (2016a); Garcia et al. (2016).

A similar transport mechanism is believed to act in the F-layer ionosphere. Here depletions in the plasma density or “bubbles” are observed in night-side equatorial regions. The rising plasma depletions are thought to trigger turbulent flows in otherwise stable regions and lead to the equatorial spread-F phenomenon, which may significantly affect the performance and reliability of radio frequency transmissions Cohen and Bowles (1961); Woodman and La Hoz (1976); Ott (1978); Hysell and Burcham (1998); Hysell (2000); Hysell and Shume (2002); Woodman (2009). Measurements of plasma depletions have also been reported from magnetically confined plasmas although their contribution to transport of plasma is still debated Boedo et al. (2003); Cheng et al. (2010); Nold et al. (2010).

In this contribution scrape-off layer plasmas as well as ionospheric plasmas are modeled by drift-fluid equations where we ignore magnetic field inhomogenity for the latter one. This simplification results in incompressible flows. As noted in Kube et al. (2016b), compressible drifts significantly alter the dynamics of seeded perturbations with low peak amplitudes relative to the background level. We further discuss the effect of the seeded perturbations’ inertial mass on the acceleration of the structure Kendl (2015). Using the conservation laws of the model equations we derive an expression that relates the acceleration of pressure perturbations to its initial amplitude relative to the background. An empirical model is proposed that is shown to reproduce velocities and accelerations taken from numerical simulations over a broad range of initial density amplitudes.

In drift-fluid models the continuity equation

[TABLE]

describes the dynamics of the electron density $n$. Here ${\mathbf{u}}_{E}:=(\hat{{\mathbf{b}}}\times\nabla\phi)/B$ gives the electric drift velocity in a magnetic field ${\mathbf{B}}:=B\hat{{\mathbf{b}}}$ and an electric potential $\phi$. We neglect contributions of the diamagnetic drift Kube et al. (2016b).

Equation (1) is closed by invoking quasineutrality, i.e. the divergence of the ion polarization, the electron diamagnetic and the gravitational drift currents must vanish

[TABLE]

Here we denote $\nabla_{\perp}\phi/B:=-\hat{{\mathbf{b}}}\times{\mathbf{u}}_{E}$, the electron diamagnetic drift ${\mathbf{u}}_{d}:=-T_{e}(\hat{{\mathbf{b}}}\times\nabla n)/enB$ with the electron temperature $T_{e}$, the ion gravitational drift velocity ${\mathbf{u}}_{g}:=m_{i}\hat{{\mathbf{b}}}\times{\mathbf{g}}/B$ with ion mass $m_{i}$, and the ion gyro-frequency $\Omega:=eB/m_{i}$.

Combining Eq. (2) with Eq. (1) yields

[TABLE]

with the polarization charge density $\rho=\nabla\cdot(n\nabla_{\perp}\phi/\Omega B)$ and ${\mathbf{u}}_{\psi}:=\hat{{\mathbf{b}}}\times\nabla\psi/B$ with $\psi:=m_{i}{\mathbf{u}}_{E}^{2}/2e$. We exploit this form of Eq. (2) in our numerical simulations.

Equations (1) and (2) respectively (3) have several invariants. First, in Eq. (1) the relative particle number $M(t):=\int\mathrm{dA}\,(n-n_{0})$ is conserved over time $\mathrm{d}M(t)/\mathrm{d}t=0$. Furthermore, we integrate $(T_{e}(1+\ln n)-T_{e}\ln B)\partial_{t}n$ as well as $-e\phi\partial_{t}\rho-(m_{i}{\mathbf{u}}_{E}^{2}/2+gm_{i}x-T_{e}\ln B)\partial_{t}n$ over the domain to get, disregarding boundary contributions,

[TABLE]

where we define the entropy $S(t):=\int\mathrm{dA}\,[n\ln(n/n_{0})-(n-n_{0})]$, the kinetic energy $E(t):=m_{i}\int\mathrm{dA}\,n{\mathbf{u}}_{E}^{2}/2$ and the potential energies $G(t):=m_{i}g\int\mathrm{dA}\,x(n-n_{0})$ and $H(t):=T_{e}\int\mathrm{dA}\,(n-n_{0})\ln(B^{-1})$. Note that $n\ln(n/n_{0})-n+n_{0}\approx(n-n_{0})^{2}/2$ for $|(n-n_{0})/n_{0}|\ll 1$ and $S(t)$ thus reduces to the local entropy form in Reference Kube et al. (2016b).

We now set up a gravitational field ${\mathbf{g}}=g\hat{x}$ and a constant homogeneous background magnetic field ${\mathbf{B}}=B_{0}\hat{z}$ in a Cartesian coordinate system. Then the divergences of the electric and gravitational drift velocities $\nabla\cdot{\mathbf{u}}_{E}$ and $\nabla\cdot{\mathbf{u}}_{g}$ and the diamagnetic current $\nabla\cdot(n{\mathbf{u}}_{d})$ vanish, which makes the flow incompressible. Furthermore, the magnetic potential energy vanishes $H(t)=0$.

In a second system we model the inhomogeneous magnetic field present in tokamaks as ${\mathbf{B}}:=B_{0}(1+x/R_{0})^{-1}\hat{z}$ and neglect the gravitational drift ${\mathbf{u}}_{g}=0$. Then, the potential energy $G(t)=0$. Note that $H(t)=m_{i}C_{\mathrm{s}}^{2}/R_{0}\int\mathrm{dA}\,x(n-n_{0})+\mathcal{O}(R_{0}^{-2})$ reduces to $G(t)$ with the effective gravity $g_{\text{eff}}:=C_{\mathrm{s}}^{2}/R_{0}$ with $C_{\mathrm{s}}^{2}:=T_{e}/m_{i}$. For the rest of this letter we treat $g$ and $g_{\text{eff}}$ as well as $G(t)$ and $H(t)$ on the same footing. The magnetic field inhomogeneity thus entails compressible flows, which is the only difference to the model describing dynamics in a homogeneous magnetic field introduced above. Since both $S(t)\geq 0$ and $E(t)\geq 0$ we further derive from Eq. (4) and Eq. (5) that the kinetic energy is bounded by $E(t)\leq T_{e}S(t)+E(t)=T_{e}S(0)$; a feature absent from the gravitational system with incompressible flows, where $S(t)=S(0)$.

We now show that the invariants Eqs. (4) and (5) present restrictions on the velocity and acceleration of plasma blobs. First, we define the blobs’ center of mass (COM) via $X(t):=\int\mathrm{dA}\,x(n-n_{0})/M$ and its COM velocity as $V(t):=\mathrm{d}X(t)/\mathrm{d}t$. The latter is proportional to the total radial particle flux Garcia et al. (2006); Held et al. (2016). We assume that $n>n_{0}$ and $(n-n_{0})^{2}/2\leq[n\ln(n/n_{0})-(n-n_{0})]n$ to show for both systems

[TABLE]

Here we use the Cauchy-Schwartz inequality and $\phi_{y}:=\partial\phi/\partial y$. Note that although we derive the inequality Eq. (6) only for amplitudes $\triangle n>0$ we assume that the results also hold for depletions. This is justified by our numerical results later in this letter. If we initialize our density field with a seeded blob of radius $\ell$ and amplitude $\triangle n$ as

[TABLE]

and $\phi({\mathbf{x}},0)=0$, we immediately have $M:=M(0)=2\pi\ell^{2}\triangle n$, $E(0)=G(0)=0$ and $S(0)=2\pi\ell^{2}f(\triangle n)$, where $f(\triangle n)$ captures the amplitude dependence of the integral for $S(0)$.

The acceleration for both incompressible and compressible flows can be estimated by assuming a linear acceleration $V=A_{0}t$ and $X=A_{0}t^{2}/2$ Held et al. (2016) and using $E(t)=G(t)=m_{i}gMX(t)$ in Eq. (6)

[TABLE]

Here, we use the Padé approximation of order $(1/1)$ of $2S(0)/M$ and define a model parameter $\mathcal{Q}$ with $0<\mathcal{Q}\leq 1$ to be determined by numerical simulations. Note that the Padé approximation is a better approximation than a simple truncated Taylor expansion especially for large relative amplitudes of order unity. Eq. (8) predicts that $A_{0}/g\sim\triangle n/n_{0}$ for small amplitudes $|\triangle n/n_{0}|<1$ and $A_{0}\sim g$ for very large amplitudes $\triangle n/n_{0}\gg 1$, which confirms the predictions in Pécseli et al. (2016) and reproduces the limits discussed in Angus and Umansky (2014).

As pointed out earlier for compressible flows Eq. (6) can be further estimated

[TABLE]

We therefore have a restriction on the maximum COM velocity for compressible flows, which is absent for incompressible flows

[TABLE]

For $|\triangle n/n_{0}|<1$ Eq. (10) reduces to the linear scaling derived in Kube et al. (2016b). Finally, a scale analysis of Eq. (3) shows that Ott (1978); Garcia et al. (2005); Held et al. (2016)

[TABLE]

This equation predicts a square root dependence of the center of mass velocity on amplitude and size.

We now propose a simple phenomenological model that captures the essential dynamics of blobs and depletions in the previously stated systems. More specifically the model reproduces the acceleration Eq. (8) with and without Boussinesq approximation, the square root scaling for the COM velocity Eq. (11) for incompressible flows as well as the relation between the square root scaling Eq. (11) and the linear scaling Eq. (10) for compressible flows. The basic idea is that the COM of blobs behaves like the one of an infinitely long plasma column immersed in an ambient plasma. The dynamics of this column reduces to the one of a two-dimensional ball. This idea is similar to the analytical “top hat” density solution for blob dynamics recently studied in Pécseli et al. (2016). The ball is subject to buoyancy as well as linear and nonlinear friction

[TABLE]

The gravity $g$ has a positive sign in the coordinate system; sgn$(f)$ is the sign function. The first term on the right hand side is the buoyancy, where $M_{\text{g}}:=\pi\ell^{2}(n_{0}+\mathcal{Q}\triangle n/2)$ is the gravitational mass of the ball with radius $\ell$ and $M_{\mathrm{p}}:=n_{0}\pi\ell^{2}$ is the mass of the displaced ambient plasma. Note that if $\triangle n<0$ the ball represents a depletion and the buoyancy term has a negative sign, i.e. the depletion will rise. We introduce an inertial mass $M_{\text{i}}:=\pi\ell^{2}(n_{0}+2\triangle n/9)$ different from the gravitational mass $M_{\text{g}}$ in order to recover the initial acceleration in Eq. (8). We interpret the parameters $\mathcal{Q}$ and $2/9$ as geometrical factors that capture the difference of the actual blob form from the idealized “top hat” solution. Also note that the Boussinesq approximation appears in the model as a neglect of inertia, $M_{\text{i}}=\pi\ell^{2}n_{0}$.

The second term is the linear friction term with coefficient $c_{1}(\ell)$, which depends on the size of the ball. If we disregard the nonlinear friction, $c_{2}=0$, Eq. (12) directly yields a maximum velocity $c_{1}V^{*}=\pi\ell^{2}ng\mathcal{Q}\triangle n/2$. From our previous considerations $\max V/C_{\mathrm{s}}=\mathcal{Q}\triangle n/2n_{0}$, we thus identify

[TABLE]

The linear friction coefficient thus depends on the gravity and the size of the ball.

The last term in (12) is the nonlinear friction. The sign of the force depends on whether the ball rises or falls in the ambient plasma. If we disregard linear friction $c_{1}=0$, we have the maximum velocity $V^{*}=\sigma(\triangle n)\sqrt{\pi\ell^{2}|\triangle n|g\mathcal{Q}/c_{2}}$, which must equal $\max V=\sigma(\triangle n)\mathcal{R}\sqrt{g\ell|\triangle n/n_{0}|}$ and thus

[TABLE]

Inserting $c_{1}$ and $c_{2}$ into Eq. (12) we can derive the maximum absolute velocity in the form

[TABLE]

and thus have a concise expression for $\max|V|$ that captures both the linear scaling (10) as well as the square root scaling (11). With Eq. (8) and Eq. (11) respectively Eq. (15) we finally arrive at an analytical expression for the time at which the maximum velocity is reached via $t_{\max V}\sim\max V/A_{0}$. Its inverse $\gamma:=t_{\max V}^{-1}$ gives the global interchange growth rate, for which an empirical expression was presented in Reference Held et al. (2016).

We use the open source library FELTOR to simulate Eqs. (1) and (3) with and without drift compression. For numerical stabilty we added small diffusive terms on the right hand sides of the equations. The discontinuous Galerkin methods employ three polynomial coefficients and a minimum of $N_{x}=N_{y}=768$ grid cells. The box size is $50\ell$ in order to mitigate influences of the finite box size on the blob dynamics. Moreover, we used the invariants in Eqs. (4) and (5) as consistency tests to verify the code and repeated simulations also in a gyrofluid model. No differences to the results presented here were found. Initial perturbations on the particle density field are given by Eq. (7), where the perturbation amplitude $\triangle n/n_{0}$ was chosen between $10^{-3}$ and $20$ for blobs and $-10^{0}$ and $-10^{-3}$ for depletions. Due to computational reasons we show results only for $\triangle n/n_{0}\leq 20$. For compressible flows we consider two different cases $\ell/R_{0}=10^{-2}$ and $\ell/R_{0}=10^{-3}$. For incompressible flows Eq. (1) and (3) can be normalized such that the blob radius is absent from the equations Ott (1978); Kube and Garcia (2012). The simulations of incompressible flows can thus be used for both sizes. The numerical code as well as input parameters and output data can be found in the supplemental dataset to this contribution Wiesenberger et al. (2017).

In Fig. 1 we plot the maximum COM velocity for blobs with and without drift compression. For incompressible flows blobs follow the square root scaling almost perfectly. Only at very large amplitudes velocities are slightly below the predicted values. For small amplitudes we observe that the compressible blobs follow a linear scaling. When the amplitudes increase there is a transition to the square root scaling at around $\triangle n/n_{0}\simeq 0.5$ for $\ell/R_{0}=10^{-2}$ and $\triangle n/n_{0}\simeq 0.05$ for $\ell/R_{0}=10^{-3}$, which is consistent with Eq. (15) and Reference Kube et al. (2016b). In the transition regions the simulated velocities are slightly larger than the predicted ones from Eq. (15). Beyond these amplitudes the velocities of compressible and incompressible blobs align.

In Fig. 2 we show the maximum radial COM velocity for depletions instead of blobs. For relative amplitudes below $|\triangle n|/n_{0}\simeq 0.5$ (right of unity in the plot) the velocities coincide with the corresponding blob velocities in Fig. 1. For amplitudes larger than $|\triangle n|/n_{0}\simeq 0.5$ the velocities follow the square root scaling. We observe that for plasma depletions beyond $90$ percent the velocities in both systems reach a constant value that is very well predicted by the square root scaling.

In Fig. 3 we show the average acceleration of blobs for compressible and incompressible flows computed by dividing the maximum velocity $\max V$ by the time to reach this velocity $t_{\max V}$. We compare the simulation results to the theoretical predictions Eq. (8) of our model with and without inertia. The results of the compressible and incompressible systems coincide and fit very well to our theoretical values. For amplitudes larger than unity the acceleration deviates significantly from the prediction with Boussinesq approximation.

In Fig. 4 we show the simulated acceleration of depletions in the compressible and the incompressible systems. We compare the simulation results to the theoretical predictions Eq. (8) of our model with and without inertia. Deviations from our theoretical prediction Eq. (8) are visible for amplitudes smaller than $\triangle n/n_{0}\simeq-0.5$ (left of unity in the plot). The relative deviations are small at around $20$ percent. As in Fig. 2 the acceleration reaches a constant values for plasma depletions of more than $90$ percent. Comparing Fig. 4 to Fig. 3 the asymmetry between blobs and depletions becomes apparent. While the acceleration of blobs is reduced for large amplitudes compared to a linear dependence the acceleration of depletions is increased. In the language of our simple buoyancy model the inertia of depletions is reduced but increased for blobs.

In conclusion we discuss the dynamics of seeded blobs and depletions in a compressible and an incompressible system. With only two fit parameters our theoretical results reproduce the numerical COM velocities and accelerations over five orders of magnitude. We derive the amplitude dependence of the acceleration of blobs and depletions from the conservation laws of our systems in Eq. (8). From the same inequality a linear regime is derived in the compressible system for ratios of amplitudes to sizes smaller than a critical value. In this regime the blob and depletion velocity depends linearly on the initial amplitude and is independent of size. The regime is absent from the system with incompressible flows. Our theoretical results are verified by numerical simulations for all amplitudes that are relevant in magnetic fusion devices. Finally, we suggest a new empirical blob model that captures the detailed dynamics of more complicated models. The Boussinesq approximation is clarified as the absence of inertia and a thus altered acceleration of blobs and depletions. The maximum blob velocity is not altered by the Boussinesq approximation.

The authors were supported with financial subvention from the Research Council of Norway under grant 240510/F20. M.W. and M.H. were supported by the Austrian Science Fund (FWF) Y398. The computational results presented have been achieved in part using the Vienna Scientific Cluster (VSC). Part of this work was performed on the Abel Cluster, owned by the University of Oslo and the Norwegian metacenter for High Performance Computing (NOTUR), and operated by the Department for Research Computing at USIT, the University of Oslo IT-department. This work has been carried out within the framework of the EUROfusion Consortium and has received funding from the Euratom research and training programme 2014-2018 under grant agreement No 633053. The views and opinions expressed herein do not necessarily reflect those of the European Commission.