  We model neutron star cooling with several microscopic nuclear equations of
state based on different nucleon-nucleon interactions and three-body forces,
and compatible with the recent GW170817 neutron star merger event. They all
feature strong direct Urca processes. We find that all models are able to
describe well the current set of cooling data for isolated neutron stars,
provided that large and extended proton 1S0 gaps and no neutron 3PF2 gaps are
active in the stellar matter. We then analyze the neutron star mass
distributions predicted by the different models and single out the preferred
ones.
