Improving SIESTA numerical convergence

SIESTA is a well known ideal MHD non linear equilibrium solver [1] capable of finding solutions that may exhibit rather complicated topologies including both magnetic islands and stochastic regions. Recently, it was extended to be capable of providing free-plasma boundary solutions as well [2]. In this work, various advances targetting the improvement of the numerical convergence of SIESTA towards the equilibrium are presented. These advances will be of use, not only for increasing the computational efficiency for cases that SIESTA already solves with ease, but to significantly improve the performance of SIESTA in configurations where convergence is slow, is plagued by plateaus of apparently no convergence or simply fail to converge. To achieve this, it has been developed a new second order advance for both p and B perturbations (currently, it is only first order) that allows to improve the equivalence between the MHD energy variation and the work associated with the force for any given displacement. By carefully discretizing the resulting advance equations, this equivalence can be made exact for the pressure term (up to machine accuracy) and second order accurate in the displacement for the magnetic field term.This equivalence has the additional associated benefit of ensuring a fully symmetric Jacobean down to machine precision. Analytical results and first numerical tests will be presented.

[1] S. P. Hirshman, R. Sanchez, and C. R. Cook, Phys. Plasmas 18, 062504 (2011)

[2] H. Peraza-Rodriguez, J. M. Reynolds-Barredo, R. Sanchez, J. Geiger, V. Tribaldos, S. P. Hirshman, M. Cianciosa, Phys. of Plasmas 24, 082516 (2017)