BPS Skyrme neutron stars in generalized gravity

We study the coupling of nuclear matter described by the BPS Skyrme model to generalized gravity. Concretely, we consider the Starobinsky model which provides the leading-order correction to the Einstein-Hilbert action. Static solutions describing neutron stars are found both for the full field theory and for the mean-field approximation. We always consider the full Starobinsky model in the nonperturbative approach, using appropriately generalized shooting methods for the numerical neutron star calculations. Many of our results are similar to previous investigations of neutron stars for the Starobinsky model using other models of nuclear matter, but there are some surprising discrepancies. The "Newtonian mass" relevant for the surface redshift, e.g., results larger than the ADM mass in our model, in contrast to other investigations. This difference is related to the particularly high stiffness of nuclear matter described by the BPS Skyrme model and offers an interesting possibility to distinguish different models of nuclear matter within generalized gravity.