Impurity transport in a mixed-collisionality stellarator plasma

TL;DR

This paper demonstrates that impurity screening, previously thought exclusive to tokamaks, can also occur in stellarator plasmas under certain collisionality conditions, potentially improving plasma performance.

Contribution

It reveals that impurity screening in stellarators is possible in high-collisionality impurity regimes, supported by analytical and numerical methods, challenging prior assumptions.

Findings

Impurity screening can occur in stellarators with highly collisional impurities.

Analytical and numerical evidence supports the possibility of impurity outward flux.

Screening occurs in regimes where impurities are collisional and the bulk plasma is low-collisional.

Abstract

A potential threat to the performance of magnetically confined fusion plasmas is the problem of impurity accumulation, which causes the concentration of highly charged impurity ions to rise uncontrollably in the center of the plasma and spoil the energy confinement by excessive radiation. It has long been thought that the collisional transport of impurities in stellarators always leads to such accumulation (if the electric field points inwards, which is usually the case), whereas tokamaks, being axisymmetric, can benefit from "temperature screening", i.e., an outward flux of impurities driven by the temperature gradient. Here it is shown, using analytical techniques supported by results from a new numerical code, that such screening can arise in stellarator plasmas too, and indeed does so in one of the most relevant operating regimes, where the impurities are highly collisional whilst the bulk plasma is in any of the low-collisionality regimes.

Equations

$${\Gamma }_z=⟨\int f_z({\mathbf{v}}_d\cdot \nabla r)d^{3}v⟩=n_z\left(D_{11}^{zi}{A}_{1i}+D_{11}^{zz}{A}_{1z}+D_{12}^z{A}_{2i}\right),$$

$$A_{1a}=\frac{d\ln p_a}{dr}+\frac{e_a{\varphi }^{\prime }(r)}{T_a},A_{2a}=\frac{d\ln T_a}{dr}$$

$${\Gamma }_z=\frac{1}{Ze}⟨uB{R}_{zi\parallel }+(p_{z\Vert }-p_{z\perp })\frac{{\nabla }_{\Vert }(u{B}^2)}{2B}⟩.$$

$$R_{zi\parallel }=\int m_z{v}_{\parallel }{C}_{zi}(f_z,f_i)d^{3}v=m_i\int {\nu }_D^{iz}(v)v_{\parallel }{f}_i{d}^{3}v-\frac{m_i{n}_i{V}_{z\parallel }}{{\tau }_{iz}},$$

$$f^{-}\left(r,\alpha ,l,ϵ,\mu ,\sigma \right)={\int }_{l_0}^l\left[C^{+}(f)-{\mathbf{v}}_d\cdot \nabla f^{+}\right]\frac{d{l}^{\prime }}{v_{\parallel }}+X,$$

$$\overline{{\mathbf{v}}_d\cdot \nabla f^{+}}=\overline{C^{+}(f)}$$

$$\overline{{\mathbf{v}}_d\cdot \nabla \alpha }\frac{\partial F_1}{\partial \alpha }+\overline{{\mathbf{v}}_d\cdot \nabla r}\frac{\partial F_0}{\partial r}\simeq 0.$$

$$R_{zi\parallel }=\frac{p_i{m}_i}{e{\tau }_{iz}}\left(A_{1i}-\frac{3}{2}{A}_{2i}\right)(u+s)B+P(r)B-\frac{m_i{n}_i{V}_{z\parallel }}{{\tau }_{iz}},$$

$$s(l)=\frac{3}{2}{\int }_{l_{\max }}^{l}d{l}^{\prime }{\int }_{1/B_{max}}^{1/B(l^{\prime })}\frac{d\lambda }{\xi (l^{\prime })}\overline{\xi \left(\hat{\mathbf{b}}\times \nabla r\right)\cdot \nabla \left(\frac{\xi }{B}\right)}$$

$${\Gamma }_z=-\frac{m_i{p}_i}{Z{e}^2{\tau }_{iz}}[\frac{A_{1z}}{Z}\left(⟨u^2{B}^2⟩-\frac{{⟨u{B}^2⟩}^2}{⟨B^2⟩}\right)$$

$$-\left(A_{1i}-\frac{3}{2}{A}_{2i}\right)\left(⟨u\left(u+s\right)B^2⟩-⟨\left(u+s\right)B^2⟩\frac{⟨u{B}^2⟩}{⟨B^2⟩}\right)],$$

$$D_{\mathrm{PS}}=\frac{m_i{T}_i}{e^2{\tau }_{iz}}\left[⟨u^2{B}^2⟩-\frac{{⟨u{B}^2⟩}^2}{⟨B^2⟩}\right]=\frac{{\rho }_i^2}{{\tau }_{iz}}\frac{⟨J_{\Vert }^2⟩⟨B^2⟩-{⟨J_{\Vert }B⟩}^2}{(dp/dr{)}^2},$$

Full text

The transport of impurity ions in a tokamak or stellarator plasma is governed by turbulent and “neoclassical” processes, the latter being caused by the random walk executed by these particles as they travel along complicated orbits set by the magnetic-field geometry whilst colliding with other particles Hirshman ; helandersigmar2002 ; helander2014 . The turbulent transport often dominates, but the neoclassical transport can be very significant for heavy impurities in both types of device. Moreover, it is usually in the direction of the bulk-ion density gradient, i.e. inward, into the core of the plasma. To make things worse, this inward neoclassical transport has been predicted to be particularly strong, practically inevitable, in stellarators. This has been an issue of great concern for several decades. The aim of the present Letter is to show that the situation is less serious than previously thought.

The neoclassical impurity flux is of the form Hirshman ; Igitkhanov

[TABLE]

where $z$ refers to impurity ions, $i$ to the bulk (hydrogenic) ions, ${\bf v}_{d}$ denotes the drift velocity, $r$ is an arbitrary label of the magnetic surfaces serving as radial coordinate, angular brackets indicate an average over such surfaces, and

[TABLE]

denote the “thermodynamic forces”. Here $p_{a}=n_{a}T_{a}$ is the pressure of species $a$, $e_{a}$ its charge, and $\phi(r)$ the electrostatic potential, which like the density $n_{a}$ and temperature $T_{a}$ is approximately constant on magnetic surfaces. Different ion species usually have the same temperature, so for simplicity we have $A_{2i}=A_{2z}$. In axisymmetric magnetic fields, the largest transport coefficients are $D_{11}^{zi}$ and $D_{12}^{z}$, where the former is always positive and so tends to drive impurities into the plasma if $dp_{i}/dr<0$. The coefficient $D_{12}^{z}$ is usually negative and can (depending on the collision frequency) exceed $D_{11}^{zi}$ Connor . If the ion temperature profile is sufficiently steep, outward impurity transport will then result. This beneficial property is referred to as temperature screening. In tokamaks, the transport coefficients also have the property that the sum of all terms containing the radial electric field, ${\bf E}=-\phi^{\prime}(r)\nabla r$ vanish.

However, for heavy impurities ($a=z$, $e_{a}=Ze$, $m_{z}/m_{i}=O(Z)\gg 1$) in a stellarator the picture has been pessimistic, because the radial electric field term does contribute to the transport. Since this term is multiplied by the large number $Z\gg 1$ in $A_{1z}$, it tends to dominate and will drive the impurities in the direction of the electric field, which typically points inward in fusion-relevant high-density plasmas. Unless there is very strong turbulence, which impairs energy confinement, any high-$Z$ impurity will thus accumulate in the core, as is indeed often observed in experiments W7-A ; Igitkhanov ; Hirsch . A notable exception occurs in low-density “impurity-hole” plasmas in the Large Helical Device Ida .

These conclusions were mostly based on calculations that approximate the collisions between different particles species by a simple scattering operator, sometimes augmented by a term ensuring momentum conservation. This treatment is adequate when all ions are at low collisionality (set by the ratio of the major radius to the mean free path and denoted by $\nu_{*}^{ab}$ for collisions between two species $a$ and $b$). However, the regime of greatest practical interest is the one where the mean free path of the bulk ions is long but that of high-$Z$ impurities relatively short. (It scales as $Z^{-4}$ at constant impurity density.) This mixed-collisionality regime cannot be treated correctly by numerical codes that neglect interaction other than collisional scattering between different ion species. For several decades, such codes have been the work horses for neoclassical transport calculations in stellarators beidleretal2011 .

We therefore consider the neoclassical transport in a mixed-collisionality plasma with a single, highly charged impurity species, and calculate the cross-field flux $\Gamma_{z}$ by solving the kinetic equation for the distribution function $f_{z}$ both analytically and numerically, using the full Landau collision operator for the impurities. The numerical calculation can be done for arbitrary collisionality, but the analytical treatment is only possible by considering the asymptotic limit of short and long mean free paths, respectively, for the impurities and bulk ions. More details of both calculations will be published separately.

Analytical calculation. In order to evaluate the impurity particle flux analytically, it is useful to decompose it into a sum of contributions sugamanishimura2002 ; braunhelander2010 , driven by the friction against the background bulk ions and the pressure anisotropy,

[TABLE]

Here, the second term on the right is relatively small for a highly collisional population, with $T_{z}=T_{i}\equiv T$, and can be neglected when the collisionalities satisfy $\nu_{*}^{iz}\nu_{*}^{zz}\gg n_{z}Z^{1/2}/n_{i}$, a condition we take to hold in the analytic calculation. The remaining term contains the function $u$, which is related to the parallel (to the magnetic field) plasma current and pressure by $u=J_{\|}/(Bp^{\prime}(r))$, and the friction force

[TABLE]

where the deflection frequency is $\nu_{D}^{iz}(v)=(3\pi^{1/2}/4\tau_{iz})(v_{Ti}/v)^{3}$ and the collision time $\tau_{iz}=3(2\pi)^{3/2}\sqrt{m_{i}}T^{3/2}\epsilon_{0}^{2}/(n_{z}Z^{2}e^{4}\ln\Lambda)$ helandersigmar2002 . To evaluate the first term appearing in Eq. (3), we require the solution of the kinetic equation for the distribution function $f_{i}$, specifically the piece which is odd in the parallel velocity, $v_{\parallel}=\sigma|v_{\parallel}|$.

A recently developed formulation allows a unified treatment of this problem throughout the low-collisionality regimes of the bulk ions helanderetalsub . If the distribution function is split into even and odd pieces, $f^{\pm}$, the kinetic equation becomes $v_{\parallel}\nabla_{\parallel}f^{\mp}=C^{\pm}(f)-{\bf v}_{d}\cdot\nabla f^{\pm}$, where the independent coordinates are $\left(r,\alpha,l,\epsilon,\mu,\sigma\right)$. Here $\alpha$ is a label for the different field lines on the same flux surface, $l$ the arc length along the magnetic field, $\epsilon=m_{i}v^{2}/2+e_{i}\phi({\bf x})$ the energy, and $\mu=m_{i}v_{\perp}^{2}/2B$ the magnetic moment. The electrostatic potential can be set to $\phi=0$ on the surface of interest, and particle orbits are defined as passing (able to move over a whole flux surface) or magnetically trapped depending on whether the parameter $\lambda=\mu/\epsilon$ is less or greater than $1/B_{\rm max}(r)$, where $B_{\rm max}(r)$ is the maximum field strength on the surface. The odd part of the bulk-ion distribution function can thus be written as

[TABLE]

where $X\left(r,\alpha,\epsilon,\mu,\sigma\right)$ denotes an integration constant that is independent of $\alpha$ in the passing region and vanishes in the trapped region if $l_{0}$ is chosen to be a bounce point. At such points $\lambda=1/B$, the parallel velocity vanishes, and so therefore must $f^{-}$.

The parallel streaming term in the kinetic equation is annihilated by the orbit average, which is defined as a time average taken along a trajectory following the magnetic field between two consecutive bounce points for trapped particles, or many times around the torus for passing ones. (Formally, this average is obtained by multiplying the kinetic equation by $dl/v_{\|}$ and integrating.) This gives the equation

[TABLE]

for $f^{+}$, where the orbit average is indicated by $\overline{(\cdots)}$. For passing particles and for trapped ones at moderately low bulk-ion collisionality (the $1/\nu$-regime), the right-hand side dominates, and the solution is approximately Maxwellian, $f^{+}\simeq F_{0}(r,\epsilon)$. At lower bulk-ion collisionality (the $\sqrt{\nu}$-regime), trapped-particle drifts produce strong deviations from such an equilibrium and prevent effective plasma confinement unless the magnetic field is optimized to be nearly omnigeneous calvoetal2016 , or the radial electric field is sufficiently strong to produce an in-surface ${\bf E}\times{\bf B}$ drift which averages out the magnetic drift motion over an orbit hokulsrud1987 . In both cases, particles stay close to a flux surface on an orbit average, so that $f^{+}=F_{0}+F_{1}$ is determined by

[TABLE]

This equation holds for all trapped particles except those in a thin layer around the trapped-passing boundary, which is unimportant for our present purposes but regulates the transport of the bulk ions hokulsrud1987 . We thus conclude that ${\bf v}_{d}\cdot\nabla f^{+}\simeq({\bf v}_{d}\cdot\nabla r)\partial F_{0}/\partial r$ in the $1/\nu$-regime and ${\bf v}_{d}\cdot\nabla f^{+}\simeq({\bf v}_{d}\cdot\nabla r-\overline{{\bf v}_{d}\cdot\nabla r})\partial F_{0}/\partial r$ in the $\sqrt{\nu}$-regime helanderetalsub .

We now evaluate the friction force between the main ions and the impurities by substituting the solution (4) in Eq. (3). The collision term in Eq. (4) is small in the $\sqrt{\nu}$-regime but generally produces a contribution to the friction in the $1/\nu$-regime. However, a commonly used, and successful landremanetal2014 , approximation to the collision operator for transport applications takes the form of a pitch-angle scattering operator with a momentum-conserving term, in which case this contribution vanishes. We will therefore neglect it, but its impact can be tested in the final numerical comparison. The remaining terms give

[TABLE]

where terms that vary as $B$ over the flux surface have been combined into $P(r)$, a flux surface function that turns out not to affect the impurity flux. In deriving this result, we have written ${\bf v}_{d}\cdot\nabla r=\Omega_{i}^{-1}v_{\parallel}\left({\bf B}\times\nabla r\right)\cdot\nabla\left(v_{\parallel}/B\right)$, and defined $s=0$ in the $1/\nu$-regime and

[TABLE]

in the $\sqrt{\nu}$-regime, with $\xi=\sqrt{1-\lambda B}$.

The impurity flow $V_{z\parallel}$ takes the well-known form $n_{z}V_{z\parallel}=(p_{z}/Ze)A_{1z}uB+K_{z}(r)B$, where $K_{z}$ is an integration constant determined by force balance along the magnetic field, which reduces to $\left\langle BR_{zi\parallel}\right\rangle=0$ in the collisional limit hirshmansigmar1977 . Using this relation to determine $K_{z}$, we eliminate the flux function $P(r)$ from the required flux surface average $\left\langle uBR_{zi\parallel}\right\rangle$ and obtain the final expression for the impurity flux (2),

[TABLE]

where we recognize the Pfirsch-Schlüter diffusion coefficient

[TABLE]

with $\rho_{i}^{2}=m_{i}T_{i}/e^{2}\left\langle B^{2}\right\rangle$.

The impurity transport coefficients $D_{ij}^{za}$ appearing in Eq. (1) can now be identified, and we note the following points. The Schwartz inequality implies $D_{\rm PS}\geq 0$, so $D_{11}^{zz}$ is negative, as required by entropy considerations. In the $1/\nu$-regime (moderate bulk ion collisionality) $s=0$, and we see $D_{11}^{zi}=-ZD_{11}^{zz}$, as is also true in the very-high-collisionality limit where both ion species are collisional braunhelander2010 . The direct drive of the flux by the electric field thus cancels out in both of these regimes. Furthermore, we note that $D_{12}^{z}=-(3/2)D_{11}^{zi}$, so there is temperature screening in the $1/\nu$-regime. At lower bulk ion collisionality, the electric-field drive no longer cancels exactly, becoming proportional to the quantity $\left\langle usB^{2}\right\rangle-\left\langle sB^{2}\right\rangle\left\langle uB^{2}\right\rangle/\left\langle B^{2}\right\rangle$, which must be evaluated numerically. The same additional contribution, due to the trapped particle drift, appears in the transport coefficient multiplying $dT/dr$, but the relation $D_{12}^{z}=-(3/2)D_{11}^{zi}$ continues to hold, so there is either temperature screening or outward transport due to the bulk-ion density gradient, depending on the sign of the last term in Eq. (6).

This analytical calculation thus shows that (i) there can be temperature screening both in the $1/\nu$- and $\sqrt{\nu}$-regimes and (ii) the radial electric field may only weakly drive impurity transport. Both results are at odds with conventional wisdom.

Numerical calculation. A novel computational tool has recently been developed, the continuum $\delta f$ SFINCS code, which solves the coupled first-order drift-kinetic equations in general magnetic geometries and calculates the neoclassical transport for an arbitrary number of species, retaining the full linearized multi-species Landau collision operator. The numerical implementation, which includes the calculation of perturbed Rosenbluth potentials, is detailed in landremanetal2014 . SFINCS was used by Mollén et. al. mollenetal2015 in an extensive study of the transport of impurities and their effect on the bootstrap current in the recently completed Wendelstein 7-X stellarator klingeretal2017 . Puzzling indications of temperature screening were seen already in this work, which we are now able to understand and clarify in terms of the analytical theory above and further numerical results.

Figure 1 shows an example where the transport coefficients for commonly occurring C6+ impurity ions in a bulk H+ plasma were studied over a wide range of impurity collisionality, defined as $\nu_{\ast}^{zz}=R/(v_{Tz}\tau_{zz})$ with $\tau_{zz}=(m_{z}/m_{i})^{1/2}Z^{-2}\tau_{iz}$ and $R=(G+\iota I)/B_{00}$, where the magnetic field is ${\bf B}=K\nabla\psi+I\nabla\theta+G\nabla\varphi$ in Boozer coordinates and $B_{00}$ denotes the $m=n=0$ Fourier harmonic of the field strength $B(\psi,\theta,\varphi)$. In both Figs. 1a) and 1b), $E_{r}=0$ and, as expected from the analytical calculation, the coefficient $D_{11}^{zz}+D_{11}^{zi}+D_{12}^{z}$ multiplying the temperature gradient is negative, implying temperature screening up to the collisionality where the bulk ions leave the long-mean-free-path regime and become collisional, at which point temperature screening disappears braunhelander2010 . If the collision operator is replaced by pure pitch-angle scattering, however, the temperature screening disappears from the entire range of high impurity collisionality. The range of collisional impurities ($\nu_{*}^{zz}>1$) and collisionless ($\nu_{*}^{iz}<1$) bulk ions is demarked by the thick black line, which gives the value of the coefficient obtained from Eq. (6). Note that this value matches the results from both SFINCS (with the full Landau collision operator) and the DKES code if a momentum-conserving term is added to the scattering operator maassberg . However, the latter treatment fails in the very-high-collisionality limit where all species are collisional. This limit can be treated analytically braunhelander2010 and the asymptote is shown as a dashed red line, along with that for pitch-angle scattering without momentum correction, matching the corresponding numerical results in the appropriate limit.

The coefficient $D_{11}^{zz}+D_{11}^{zi}+D_{12}^{z}$ is negative also at low impurity collisionalities, but its contribution to the transport is expected to be overwhelmed by the inward flux produced by a radial electric field. The corresponding electric-field transport coefficient (the field is assumed to be weak enough to allow the $1/\nu$-regime over the collisionality range) is shown in Fig. 1b) and can be seen to be large at low collisionality but small in the mixed collisionality regime, enabling temperature screening.

Figure 2 shows transport coefficients similar to Fig. 1 but for the heavy trace impurity Fe16+ and $Z_{\rm eff}=1.07$. Again, the analytical prediction (6) is confirmed, but now there is temperature screening also in the very-high-collisionality regime where both ion species are collisional, as expected for very clean plasmas Rutherford .

When a larger radial electric field is present, the terms in Eq. (6) involving the quantity $s$ are important and can be expected to reduce the transport driven by the bulk-ion density and temperature gradients helanderetalsub . Temperature screening is still possible, but the contributions from the radial electric field in these terms no longer exactly cancel that from the impurity thermodynamic force $A_{1z}$. On the one hand, there is thus a reduction of the overall neoclassical transport (making it easier, for instance, for turbulence to expel the impurities), but on the other hand some direct electric-field-driven transport now remains. The net transport thus depends sensitively on the collisionality, charge number and the relative size of the density and temperature gradients. By way of example, Fig. 3 shows the transport coefficients in the same magnetic W7-X configuration as for Figs. 1 and 2, but in the presence of a finite radial electric field of $E_{r}=-5\rm kV/m$. Temperature screening still occurs (if the temperature gradient is strong enough), but at a lower level, as one would expect from the theory above, making the coefficient comparable to the electric-field coefficient at high collisionality.

Conclusions. In summary, we have found that highly charged, collisional impurity ions in stellarators can experience neoclassical temperature screening if the mean free path of the bulk ions is long and the temperature profile is sufficiently steep, just as in a tokamak. Under these conditions, which are common in experiments, the radial electric field must compete to drive the impurities inward. Impurity accumulation is therefore not inevitable, and impurities may be expected to enjoy outward collisional transport when the conditions are right.