kobra: a new Vlasov code for plasma wall modeling

kobra: a new Vlasov code for plasma wall modeling APPLIED PHYSICS –NUCLEAR FUSION Sebastian Konewko, Nathan Maestracci, Prof. Sven Van Loo Benchmarks: Background For the operation of tokamaks, it is fundamental to understand the interaction of the plasma with the reactor wall.In this region, i.e. the scrape off layer (SOL), particle distributions are often nonMaxwellian and particle density is sparse. Particle-in-cell codes offer astraight-forward formulation to addressing kinetic problems but in the SOL, this approach leads to statistical errors. Our goal is to develop aVlasov code, kobra,devoid of these errors and applicable to plasma wall interactions. Contact [email protected] nuclearfusion.ugent.be Sources [1] S. A. E. G. Falle. “A numerical scheme for multifluid magnetohydrodynamics”. In: Monthly Notices of the Royal Astronomical Society 344.4 (Oct. 2003), DOI: 10.1046/j.1365-8711.2003.06908.x. [2] J.A.F. Hittinger and J.W. Banks. “Block-structured adaptive mesh refinement algorithms for Vlasov simulation”. In: Journal of Computational Physics 241 (2013) DOI: https://doi.org/10.1016/j.jcp.2013.01.030. [3] Van Loo, Sven and Falle, S. A. E. G. and Hartquist, T. W. “Generation of density inhomogeneities by magnetohydrodynamic waves in two dimensions” [4] Vogman, Genia, et al. “Fourth-Order Conservative Vlasov-Maxwell Solver for Cartesian and Cylindrical Phase Space Coordinates.”, University of California, 2016. [5] R. Chodura. “Plasma–wall transition in an oblique magnetic field”. In: The Physics of Fluids 25.9 (Sept. 1982), pp. 1628–1633. ISSN: 0031-9171. DOI: 10.1063/1.863955. Numerical discretization Next Steps: •Add oblique external magnetic field to plasma sheath model •Promote to fully 4th order code •Test higher dimensional problems 1d3v, 2d2v… •Extend to Vlasov-Maxwell system Conclusions: Plasma Wall interaction: 𝜕𝑓! 𝜕t +v 𝜕𝑓! 𝜕x +𝑞! 𝑚!v×𝐵+𝐸 𝜕𝑓! 𝜕v = 0 −𝜕"𝜙 𝜕x"= 0 # n#q#E = −∇𝜙 kobra is a second order code which solves the Vlasov-Poisson system. We use an adaptive mesh routine for computational gain. Vlasov Equation (hyperbolic solver ): •Runge-Kutta midpoint scheme for temporal discretization •Fluxes are calculated with an upwind scheme with a flux limiter [1] Poisson Equation (elliptic solver): •We generate a second auxiliary grid in the spatial coordinates to calculate the scalar potential •We use a multigrid approach with a Jacobi relaxation scheme To verify that our code is producing physically meaningful results, we benchmark our progress with the strong landau damping and two stream instability. These problems have semi-analytic solutions for the growth and damping rates which we can use to validate our results. 𝑓 !=𝑣" 2𝜋 1 + 𝛼 cos 𝑥 2exp −𝑣" 2 𝑓 !=1 2𝜋 1 − 𝛼 cos 𝑥 2exp −𝑣" 2 0 5 10 15 20 25 30 Time 10°8 10°6 10°4 10°2 Electric Energy Growth Rate: 0.50 Uniform Grid Adaptive Grid 0 20 40 60 80 100 Time °1 0 1 2 3 4 Percent Deviation Uniform Grid Adaptive Grid 4th Order 0 10 20 30 40 50 Time 10°7 10°5 10°3 10°1 Electric Energy Damping Rate: °0.5904 Growth Rate: 0.1688 Uniform Grid Adaptive Grid Adaptive Mesh Refinement (AMR): •We refine cells on a cell-by-cell basis •We use the first two truncation error terms as our refinement criteria [2] following the architecture of [3] •We employ defragmentation and load balancing routines to streamline the computational efficiency Landau Damping: Two Stream Instability: 𝛼 = 0.5, 𝑥, 𝑣 ∈ 0, 4𝜋 x −10,10 𝛼 = 0.01, 𝑥, 𝑣 ∈ 0, 4𝜋 x −10 ,10 Above, a snapshot of the Landau damping and two stream instability taken at t=20. It reproduces the filamentation in the Landau damping and the particle vortex for the two stream instability. Simulations were run with 256 x 256 resolution. On the right, a simulation of the two stream instability using 7 levels of AMR with an effective resolution of 512 x 512. The damping and growth rates of both the Landau damping and two stream instability test cases: we compare both uniform and adaptive grids both with an effective resolution of 256 x 256. Analytic growth rates taken from [4]. Simulations run with AMR take less than half of the uniform grid runtime. On the right, a plot showing the total energy of the two stream system. The second-order kobra struggles with energy conservation. We are working towards a fourth order version of the code and preliminary results already seem promising! To further test and build the capabilities of kobra, we model a plasma sheath [5]. We employ an inflow boundary on the left and an absorbing wall on the right. For the Poisson equation, we use a floating-point boundary potential on each side. We begin with an empty domain, assume cold ions, and a mass ratio of 𝑚!=1836 𝑚". Our inflow conditions are: 𝑓 #! =𝑛#! 2𝜋 exp −𝑣" 2 𝑓 $! = 𝑛%!𝛿 𝑣 − 𝑣! Above, a converged plasma sheath (t=1000) with 𝑣#= 0.2 and a normalization ratio 𝑛!# = 1.115 𝑛"#. The resolution is 𝑥, 𝑣 = 128 384 . Above is the converged boundary charge plotted over time. An unchanging zero charge at the inflow side is confirmation that our system is charge neutral at the bulk plasma. We show the electron profile taken at x = 0 and x = 7. We reproduce the analytic cutoff velocity 𝑣$%& = 2q'∅(/m'and show the profile is unchanged in the body of the sheath (x < 7) reconfirming the charge density profile. The total charge density of the simulated plasma sheath. It stays constant at zero going towards the bulk plasma then increases as ions are attracted to the negatively charged wall. We are on the right track! We reproduce the damping and growth rates of standard kinetic benchmarks and plasma sheath physics. We already demonstrate computational speeding using our AMR routine. We have a fourth order code working on a single cpu and uniform mesh. °4°2 0 2 4 6 Velocity 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 Electron Density fwall =0.76 vcutoff =°1.24 ne0 p2pexp≥°v2 2¥ fe(x=0,v) fe(x=7,v) 0 100 200 300 400 500 Time °1.0 °0.8 °0.6 °0.4 °0.2 0.0 Boundary Charge Inflow Charge (left) Wall Charge (right) 0 2 4 6 8 10 X 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 Charge Density ne+ni