  We propose a high-order adaptive numerical solver for the semilinear elliptic
boundary value problem modelling magnetic plasma equilibrium in axisymmetric
confinement devices. In the fixed boundary case, the equation is posed on
curved domains with piecewise smooth curved boundaries that may present
corners. The solution method we present is based on the hybridizable
discontinuous Galerkin method and sidesteps the need for geometry-conforming
triangulations thanks to a transfer technique that allows to approximate the
solution using only a polygonal subset as computational domain. Moreover, the
solver features automatic mesh refinement driven by a residual-based a
posteriori error estimator. As the mesh is locally refined, the computational
domain is automatically updated in order to always maintain the distance
between the actual boundary and the computational boundary of the order of the
local mesh diameter. Numerical evidence is presented of the suitability of the
estimator as an approximate error measure for physically relevant equilibria
with pressure pedestals, internal transport barriers, and current holes on
realistic geometries.
