Effect of time step, neutral–neutral collisions, and an underrelaxation scheme on the numerical convergence of SOLPS-ITER plasma boundary simulations with kinetic neutrals

The numerical convergence of plasma boundary simulations with kinetic neutrals in SOLPS-ITER is studied by performing a cross-parameter scan on the time step Δt and the number of Monte Carlo particles per iteration P. A low-power D-only DEMO case, a high-power D+He+Ar DEMO case and a D+He artificial slab case are used to perform the analysis. The numerical scans reveal that for a given time step, the bias error scales with 1 / P, as expected from theory. For a given P, the bias error scales with Δt1/n, with n ∈ [ 1 , 2 ]. It is found that the D+He+Ar DEMO case requires a large amount of impurity Monte Carlo particles per iteration to achieve acceptable numerical errors. Furthermore, it is explained that neutral–neutral collisions introduce an additional bias contribution, but that it is negligible compared to the bias due to the plasma neutral coupling for the studied cases. Additionally, it is shown that global particle imbalances are linked to the numerical convergence. Finally, it is shown that an underrelaxation scheme for the kinetic neutral source terms is able to reduce the bias errors by one to two orders of magnitude.