  Quiescent hard X-ray and soft gamma-ray emission from neutron stars
constitute a promising frontier to explore axion-like-particles (ALPs). ALP
production in the core peaks at energies of a few keV to a few hundreds of keV;
subsequently, the ALPs escape and convert to photons in the magnetosphere. The
emissivity goes as $\sim T^6$ while the conversion probability is enhanced for
large magnetic fields, making magnetars, with their high core temperatures and
strong magnetic fields, ideal targets for probing ALPs. We compute the energy
spectrum of photons resulting from conversion of ALPs in the magnetosphere and
then compare it against hard X-ray data from NuSTAR, INTEGRAL, and XMM-Newton,
for a set of eight magnetars for which such data exists. Upper limits are
placed on the product of the ALP-nucleon and ALP-photon couplings. For the
production in the core, we perform a calculation of the ALP emissivity in
degenerate nuclear matter modeled by a relativistic mean field theory. The
reduction of the emissivity due to improvements to the one-pion exchange
approximation is incorporated, as is the suppression of the emissivity due to
proton superfluidity in the neutron star core. A range of core temperatures is
considered, corresponding to different models of the steady heat transfer from
the core to the stellar surface. For the subsequent conversion, we solve the
coupled differential equations mixing ALPs and photons in the magnetosphere.
The conversion occurs due to a competition between the dipolar magnetic field
and the photon refractive index induced by the external magnetic field.
Semi-analytic expressions are provided alongside the full numerical results. We
also present an analysis of the uncertainty on the axion limits we derive due
to the uncertainties in the magnetar masses, nuclear matter equation of state,
and the proton superfluid critical temperature.
